1 /*M/////////////////////////////////////////////////////////////////////////////////////// 2 // 3 // IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING. 4 // 5 // By downloading, copying, installing or using the software you agree to this license. 6 // If you do not agree to this license, do not download, install, 7 // copy or use the software. 8 // 9 // 10 // License Agreement 11 // For Open Source Computer Vision Library 12 // 13 // Copyright (C) 2000-2008, Intel Corporation, all rights reserved. 14 // Copyright (C) 2009, Willow Garage Inc., all rights reserved. 15 // Third party copyrights are property of their respective owners. 16 // 17 // Redistribution and use in source and binary forms, with or without modification, 18 // are permitted provided that the following conditions are met: 19 // 20 // * Redistribution's of source code must retain the above copyright notice, 21 // this list of conditions and the following disclaimer. 22 // 23 // * Redistribution's in binary form must reproduce the above copyright notice, 24 // this list of conditions and the following disclaimer in the documentation 25 // and/or other materials provided with the distribution. 26 // 27 // * The name of the copyright holders may not be used to endorse or promote products 28 // derived from this software without specific prior written permission. 29 // 30 // This software is provided by the copyright holders and contributors "as is" and 31 // any express or implied warranties, including, but not limited to, the implied 32 // warranties of merchantability and fitness for a particular purpose are disclaimed. 33 // In no event shall the Intel Corporation or contributors be liable for any direct, 34 // indirect, incidental, special, exemplary, or consequential damages 35 // (including, but not limited to, procurement of substitute goods or services; 36 // loss of use, data, or profits; or business interruption) however caused 37 // and on any theory of liability, whether in contract, strict liability, 38 // or tort (including negligence or otherwise) arising in any way out of 39 // the use of this software, even if advised of the possibility of such damage. 40 // 41 //M*/ 42 43 #include "precomp.hpp" 44 45 namespace cv { 46 47 Stitcher Stitcher::createDefault(bool try_use_gpu) 48 { 49 Stitcher stitcher; 50 stitcher.setRegistrationResol(0.6); 51 stitcher.setSeamEstimationResol(0.1); 52 stitcher.setCompositingResol(ORIG_RESOL); 53 stitcher.setPanoConfidenceThresh(1); 54 stitcher.setWaveCorrection(true); 55 stitcher.setWaveCorrectKind(detail::WAVE_CORRECT_HORIZ); 56 stitcher.setFeaturesMatcher(makePtr<detail::BestOf2NearestMatcher>(try_use_gpu)); 57 stitcher.setBundleAdjuster(makePtr<detail::BundleAdjusterRay>()); 58 59 #ifdef HAVE_CUDA 60 if (try_use_gpu && cuda::getCudaEnabledDeviceCount() > 0) 61 { 62 #ifdef HAVE_OPENCV_XFEATURES2D 63 stitcher.setFeaturesFinder(makePtr<detail::SurfFeaturesFinderGpu>()); 64 #else 65 stitcher.setFeaturesFinder(makePtr<detail::OrbFeaturesFinder>()); 66 #endif 67 stitcher.setWarper(makePtr<SphericalWarperGpu>()); 68 stitcher.setSeamFinder(makePtr<detail::GraphCutSeamFinderGpu>()); 69 } 70 else 71 #endif 72 { 73 #ifdef HAVE_OPENCV_XFEATURES2D 74 stitcher.setFeaturesFinder(makePtr<detail::SurfFeaturesFinder>()); 75 #else 76 stitcher.setFeaturesFinder(makePtr<detail::OrbFeaturesFinder>()); 77 #endif 78 stitcher.setWarper(makePtr<SphericalWarper>()); 79 stitcher.setSeamFinder(makePtr<detail::GraphCutSeamFinder>(detail::GraphCutSeamFinderBase::COST_COLOR)); 80 } 81 82 stitcher.setExposureCompensator(makePtr<detail::BlocksGainCompensator>()); 83 stitcher.setBlender(makePtr<detail::MultiBandBlender>(try_use_gpu)); 84 85 return stitcher; 86 } 87 88 89 Stitcher::Status Stitcher::estimateTransform(InputArrayOfArrays images) 90 { 91 return estimateTransform(images, std::vector<std::vector<Rect> >()); 92 } 93 94 95 Stitcher::Status Stitcher::estimateTransform(InputArrayOfArrays images, const std::vector<std::vector<Rect> > &rois) 96 { 97 images.getUMatVector(imgs_); 98 rois_ = rois; 99 100 Status status; 101 102 if ((status = matchImages()) != OK) 103 return status; 104 105 if ((status = estimateCameraParams()) != OK) 106 return status; 107 108 return OK; 109 } 110 111 112 113 Stitcher::Status Stitcher::composePanorama(OutputArray pano) 114 { 115 return composePanorama(std::vector<UMat>(), pano); 116 } 117 118 119 Stitcher::Status Stitcher::composePanorama(InputArrayOfArrays images, OutputArray pano) 120 { 121 LOGLN("Warping images (auxiliary)... "); 122 123 std::vector<UMat> imgs; 124 images.getUMatVector(imgs); 125 if (!imgs.empty()) 126 { 127 CV_Assert(imgs.size() == imgs_.size()); 128 129 UMat img; 130 seam_est_imgs_.resize(imgs.size()); 131 132 for (size_t i = 0; i < imgs.size(); ++i) 133 { 134 imgs_[i] = imgs[i]; 135 resize(imgs[i], img, Size(), seam_scale_, seam_scale_); 136 seam_est_imgs_[i] = img.clone(); 137 } 138 139 std::vector<UMat> seam_est_imgs_subset; 140 std::vector<UMat> imgs_subset; 141 142 for (size_t i = 0; i < indices_.size(); ++i) 143 { 144 imgs_subset.push_back(imgs_[indices_[i]]); 145 seam_est_imgs_subset.push_back(seam_est_imgs_[indices_[i]]); 146 } 147 148 seam_est_imgs_ = seam_est_imgs_subset; 149 imgs_ = imgs_subset; 150 } 151 152 UMat pano_; 153 154 #if ENABLE_LOG 155 int64 t = getTickCount(); 156 #endif 157 158 std::vector<Point> corners(imgs_.size()); 159 std::vector<UMat> masks_warped(imgs_.size()); 160 std::vector<UMat> images_warped(imgs_.size()); 161 std::vector<Size> sizes(imgs_.size()); 162 std::vector<UMat> masks(imgs_.size()); 163 164 // Prepare image masks 165 for (size_t i = 0; i < imgs_.size(); ++i) 166 { 167 masks[i].create(seam_est_imgs_[i].size(), CV_8U); 168 masks[i].setTo(Scalar::all(255)); 169 } 170 171 // Warp images and their masks 172 Ptr<detail::RotationWarper> w = warper_->create(float(warped_image_scale_ * seam_work_aspect_)); 173 for (size_t i = 0; i < imgs_.size(); ++i) 174 { 175 Mat_<float> K; 176 cameras_[i].K().convertTo(K, CV_32F); 177 K(0,0) *= (float)seam_work_aspect_; 178 K(0,2) *= (float)seam_work_aspect_; 179 K(1,1) *= (float)seam_work_aspect_; 180 K(1,2) *= (float)seam_work_aspect_; 181 182 corners[i] = w->warp(seam_est_imgs_[i], K, cameras_[i].R, INTER_LINEAR, BORDER_CONSTANT, images_warped[i]); 183 sizes[i] = images_warped[i].size(); 184 185 w->warp(masks[i], K, cameras_[i].R, INTER_NEAREST, BORDER_CONSTANT, masks_warped[i]); 186 } 187 188 std::vector<UMat> images_warped_f(imgs_.size()); 189 for (size_t i = 0; i < imgs_.size(); ++i) 190 images_warped[i].convertTo(images_warped_f[i], CV_32F); 191 192 LOGLN("Warping images, time: " << ((getTickCount() - t) / getTickFrequency()) << " sec"); 193 194 // Find seams 195 exposure_comp_->feed(corners, images_warped, masks_warped); 196 seam_finder_->find(images_warped_f, corners, masks_warped); 197 198 // Release unused memory 199 seam_est_imgs_.clear(); 200 images_warped.clear(); 201 images_warped_f.clear(); 202 masks.clear(); 203 204 LOGLN("Compositing..."); 205 #if ENABLE_LOG 206 t = getTickCount(); 207 #endif 208 209 UMat img_warped, img_warped_s; 210 UMat dilated_mask, seam_mask, mask, mask_warped; 211 212 //double compose_seam_aspect = 1; 213 double compose_work_aspect = 1; 214 bool is_blender_prepared = false; 215 216 double compose_scale = 1; 217 bool is_compose_scale_set = false; 218 219 UMat full_img, img; 220 for (size_t img_idx = 0; img_idx < imgs_.size(); ++img_idx) 221 { 222 LOGLN("Compositing image #" << indices_[img_idx] + 1); 223 #if ENABLE_LOG 224 int64 compositing_t = getTickCount(); 225 #endif 226 227 // Read image and resize it if necessary 228 full_img = imgs_[img_idx]; 229 if (!is_compose_scale_set) 230 { 231 if (compose_resol_ > 0) 232 compose_scale = std::min(1.0, std::sqrt(compose_resol_ * 1e6 / full_img.size().area())); 233 is_compose_scale_set = true; 234 235 // Compute relative scales 236 //compose_seam_aspect = compose_scale / seam_scale_; 237 compose_work_aspect = compose_scale / work_scale_; 238 239 // Update warped image scale 240 warped_image_scale_ *= static_cast<float>(compose_work_aspect); 241 w = warper_->create((float)warped_image_scale_); 242 243 // Update corners and sizes 244 for (size_t i = 0; i < imgs_.size(); ++i) 245 { 246 // Update intrinsics 247 cameras_[i].focal *= compose_work_aspect; 248 cameras_[i].ppx *= compose_work_aspect; 249 cameras_[i].ppy *= compose_work_aspect; 250 251 // Update corner and size 252 Size sz = full_img_sizes_[i]; 253 if (std::abs(compose_scale - 1) > 1e-1) 254 { 255 sz.width = cvRound(full_img_sizes_[i].width * compose_scale); 256 sz.height = cvRound(full_img_sizes_[i].height * compose_scale); 257 } 258 259 Mat K; 260 cameras_[i].K().convertTo(K, CV_32F); 261 Rect roi = w->warpRoi(sz, K, cameras_[i].R); 262 corners[i] = roi.tl(); 263 sizes[i] = roi.size(); 264 } 265 } 266 if (std::abs(compose_scale - 1) > 1e-1) 267 { 268 #if ENABLE_LOG 269 int64 resize_t = getTickCount(); 270 #endif 271 resize(full_img, img, Size(), compose_scale, compose_scale); 272 LOGLN(" resize time: " << ((getTickCount() - resize_t) / getTickFrequency()) << " sec"); 273 } 274 else 275 img = full_img; 276 full_img.release(); 277 Size img_size = img.size(); 278 279 LOGLN(" after resize time: " << ((getTickCount() - compositing_t) / getTickFrequency()) << " sec"); 280 281 Mat K; 282 cameras_[img_idx].K().convertTo(K, CV_32F); 283 284 #if ENABLE_LOG 285 int64 pt = getTickCount(); 286 #endif 287 // Warp the current image 288 w->warp(img, K, cameras_[img_idx].R, INTER_LINEAR, BORDER_CONSTANT, img_warped); 289 LOGLN(" warp the current image: " << ((getTickCount() - pt) / getTickFrequency()) << " sec"); 290 #if ENABLE_LOG 291 pt = getTickCount(); 292 #endif 293 294 // Warp the current image mask 295 mask.create(img_size, CV_8U); 296 mask.setTo(Scalar::all(255)); 297 w->warp(mask, K, cameras_[img_idx].R, INTER_NEAREST, BORDER_CONSTANT, mask_warped); 298 LOGLN(" warp the current image mask: " << ((getTickCount() - pt) / getTickFrequency()) << " sec"); 299 #if ENABLE_LOG 300 pt = getTickCount(); 301 #endif 302 303 // Compensate exposure 304 exposure_comp_->apply((int)img_idx, corners[img_idx], img_warped, mask_warped); 305 LOGLN(" compensate exposure: " << ((getTickCount() - pt) / getTickFrequency()) << " sec"); 306 #if ENABLE_LOG 307 pt = getTickCount(); 308 #endif 309 310 img_warped.convertTo(img_warped_s, CV_16S); 311 img_warped.release(); 312 img.release(); 313 mask.release(); 314 315 // Make sure seam mask has proper size 316 dilate(masks_warped[img_idx], dilated_mask, Mat()); 317 resize(dilated_mask, seam_mask, mask_warped.size()); 318 319 bitwise_and(seam_mask, mask_warped, mask_warped); 320 321 LOGLN(" other: " << ((getTickCount() - pt) / getTickFrequency()) << " sec"); 322 #if ENABLE_LOG 323 pt = getTickCount(); 324 #endif 325 326 if (!is_blender_prepared) 327 { 328 blender_->prepare(corners, sizes); 329 is_blender_prepared = true; 330 } 331 332 LOGLN(" other2: " << ((getTickCount() - pt) / getTickFrequency()) << " sec"); 333 334 LOGLN(" feed..."); 335 #if ENABLE_LOG 336 int64 feed_t = getTickCount(); 337 #endif 338 // Blend the current image 339 blender_->feed(img_warped_s, mask_warped, corners[img_idx]); 340 LOGLN(" feed time: " << ((getTickCount() - feed_t) / getTickFrequency()) << " sec"); 341 LOGLN("Compositing ## time: " << ((getTickCount() - compositing_t) / getTickFrequency()) << " sec"); 342 } 343 344 #if ENABLE_LOG 345 int64 blend_t = getTickCount(); 346 #endif 347 UMat result, result_mask; 348 blender_->blend(result, result_mask); 349 LOGLN("blend time: " << ((getTickCount() - blend_t) / getTickFrequency()) << " sec"); 350 351 LOGLN("Compositing, time: " << ((getTickCount() - t) / getTickFrequency()) << " sec"); 352 353 // Preliminary result is in CV_16SC3 format, but all values are in [0,255] range, 354 // so convert it to avoid user confusing 355 result.convertTo(pano, CV_8U); 356 357 return OK; 358 } 359 360 361 Stitcher::Status Stitcher::stitch(InputArrayOfArrays images, OutputArray pano) 362 { 363 Status status = estimateTransform(images); 364 if (status != OK) 365 return status; 366 return composePanorama(pano); 367 } 368 369 370 Stitcher::Status Stitcher::stitch(InputArrayOfArrays images, const std::vector<std::vector<Rect> > &rois, OutputArray pano) 371 { 372 Status status = estimateTransform(images, rois); 373 if (status != OK) 374 return status; 375 return composePanorama(pano); 376 } 377 378 379 Stitcher::Status Stitcher::matchImages() 380 { 381 if ((int)imgs_.size() < 2) 382 { 383 LOGLN("Need more images"); 384 return ERR_NEED_MORE_IMGS; 385 } 386 387 work_scale_ = 1; 388 seam_work_aspect_ = 1; 389 seam_scale_ = 1; 390 bool is_work_scale_set = false; 391 bool is_seam_scale_set = false; 392 UMat full_img, img; 393 features_.resize(imgs_.size()); 394 seam_est_imgs_.resize(imgs_.size()); 395 full_img_sizes_.resize(imgs_.size()); 396 397 LOGLN("Finding features..."); 398 #if ENABLE_LOG 399 int64 t = getTickCount(); 400 #endif 401 402 for (size_t i = 0; i < imgs_.size(); ++i) 403 { 404 full_img = imgs_[i]; 405 full_img_sizes_[i] = full_img.size(); 406 407 if (registr_resol_ < 0) 408 { 409 img = full_img; 410 work_scale_ = 1; 411 is_work_scale_set = true; 412 } 413 else 414 { 415 if (!is_work_scale_set) 416 { 417 work_scale_ = std::min(1.0, std::sqrt(registr_resol_ * 1e6 / full_img.size().area())); 418 is_work_scale_set = true; 419 } 420 resize(full_img, img, Size(), work_scale_, work_scale_); 421 } 422 if (!is_seam_scale_set) 423 { 424 seam_scale_ = std::min(1.0, std::sqrt(seam_est_resol_ * 1e6 / full_img.size().area())); 425 seam_work_aspect_ = seam_scale_ / work_scale_; 426 is_seam_scale_set = true; 427 } 428 429 if (rois_.empty()) 430 (*features_finder_)(img, features_[i]); 431 else 432 { 433 std::vector<Rect> rois(rois_[i].size()); 434 for (size_t j = 0; j < rois_[i].size(); ++j) 435 { 436 Point tl(cvRound(rois_[i][j].x * work_scale_), cvRound(rois_[i][j].y * work_scale_)); 437 Point br(cvRound(rois_[i][j].br().x * work_scale_), cvRound(rois_[i][j].br().y * work_scale_)); 438 rois[j] = Rect(tl, br); 439 } 440 (*features_finder_)(img, features_[i], rois); 441 } 442 features_[i].img_idx = (int)i; 443 LOGLN("Features in image #" << i+1 << ": " << features_[i].keypoints.size()); 444 445 resize(full_img, img, Size(), seam_scale_, seam_scale_); 446 seam_est_imgs_[i] = img.clone(); 447 } 448 449 // Do it to save memory 450 features_finder_->collectGarbage(); 451 full_img.release(); 452 img.release(); 453 454 LOGLN("Finding features, time: " << ((getTickCount() - t) / getTickFrequency()) << " sec"); 455 456 LOG("Pairwise matching"); 457 #if ENABLE_LOG 458 t = getTickCount(); 459 #endif 460 (*features_matcher_)(features_, pairwise_matches_, matching_mask_); 461 features_matcher_->collectGarbage(); 462 LOGLN("Pairwise matching, time: " << ((getTickCount() - t) / getTickFrequency()) << " sec"); 463 464 // Leave only images we are sure are from the same panorama 465 indices_ = detail::leaveBiggestComponent(features_, pairwise_matches_, (float)conf_thresh_); 466 std::vector<UMat> seam_est_imgs_subset; 467 std::vector<UMat> imgs_subset; 468 std::vector<Size> full_img_sizes_subset; 469 for (size_t i = 0; i < indices_.size(); ++i) 470 { 471 imgs_subset.push_back(imgs_[indices_[i]]); 472 seam_est_imgs_subset.push_back(seam_est_imgs_[indices_[i]]); 473 full_img_sizes_subset.push_back(full_img_sizes_[indices_[i]]); 474 } 475 seam_est_imgs_ = seam_est_imgs_subset; 476 imgs_ = imgs_subset; 477 full_img_sizes_ = full_img_sizes_subset; 478 479 if ((int)imgs_.size() < 2) 480 { 481 LOGLN("Need more images"); 482 return ERR_NEED_MORE_IMGS; 483 } 484 485 return OK; 486 } 487 488 489 Stitcher::Status Stitcher::estimateCameraParams() 490 { 491 detail::HomographyBasedEstimator estimator; 492 if (!estimator(features_, pairwise_matches_, cameras_)) 493 return ERR_HOMOGRAPHY_EST_FAIL; 494 495 for (size_t i = 0; i < cameras_.size(); ++i) 496 { 497 Mat R; 498 cameras_[i].R.convertTo(R, CV_32F); 499 cameras_[i].R = R; 500 //LOGLN("Initial intrinsic parameters #" << indices_[i] + 1 << ":\n " << cameras_[i].K()); 501 } 502 503 bundle_adjuster_->setConfThresh(conf_thresh_); 504 if (!(*bundle_adjuster_)(features_, pairwise_matches_, cameras_)) 505 return ERR_CAMERA_PARAMS_ADJUST_FAIL; 506 507 // Find median focal length and use it as final image scale 508 std::vector<double> focals; 509 for (size_t i = 0; i < cameras_.size(); ++i) 510 { 511 //LOGLN("Camera #" << indices_[i] + 1 << ":\n" << cameras_[i].K()); 512 focals.push_back(cameras_[i].focal); 513 } 514 515 std::sort(focals.begin(), focals.end()); 516 if (focals.size() % 2 == 1) 517 warped_image_scale_ = static_cast<float>(focals[focals.size() / 2]); 518 else 519 warped_image_scale_ = static_cast<float>(focals[focals.size() / 2 - 1] + focals[focals.size() / 2]) * 0.5f; 520 521 if (do_wave_correct_) 522 { 523 std::vector<Mat> rmats; 524 for (size_t i = 0; i < cameras_.size(); ++i) 525 rmats.push_back(cameras_[i].R.clone()); 526 detail::waveCorrect(rmats, wave_correct_kind_); 527 for (size_t i = 0; i < cameras_.size(); ++i) 528 cameras_[i].R = rmats[i]; 529 } 530 531 return OK; 532 } 533 534 535 Ptr<Stitcher> createStitcher(bool try_use_gpu) 536 { 537 Ptr<Stitcher> stitcher = makePtr<Stitcher>(); 538 stitcher->setRegistrationResol(0.6); 539 stitcher->setSeamEstimationResol(0.1); 540 stitcher->setCompositingResol(Stitcher::ORIG_RESOL); 541 stitcher->setPanoConfidenceThresh(1); 542 stitcher->setWaveCorrection(true); 543 stitcher->setWaveCorrectKind(detail::WAVE_CORRECT_HORIZ); 544 stitcher->setFeaturesMatcher(makePtr<detail::BestOf2NearestMatcher>(try_use_gpu)); 545 stitcher->setBundleAdjuster(makePtr<detail::BundleAdjusterRay>()); 546 547 #ifdef HAVE_CUDA 548 if (try_use_gpu && cuda::getCudaEnabledDeviceCount() > 0) 549 { 550 #ifdef HAVE_OPENCV_NONFREE 551 stitcher->setFeaturesFinder(makePtr<detail::SurfFeaturesFinderGpu>()); 552 #else 553 stitcher->setFeaturesFinder(makePtr<detail::OrbFeaturesFinder>()); 554 #endif 555 stitcher->setWarper(makePtr<SphericalWarperGpu>()); 556 stitcher->setSeamFinder(makePtr<detail::GraphCutSeamFinderGpu>()); 557 } 558 else 559 #endif 560 { 561 #ifdef HAVE_OPENCV_NONFREE 562 stitcher->setFeaturesFinder(makePtr<detail::SurfFeaturesFinder>()); 563 #else 564 stitcher->setFeaturesFinder(makePtr<detail::OrbFeaturesFinder>()); 565 #endif 566 stitcher->setWarper(makePtr<SphericalWarper>()); 567 stitcher->setSeamFinder(makePtr<detail::GraphCutSeamFinder>(detail::GraphCutSeamFinderBase::COST_COLOR)); 568 } 569 570 stitcher->setExposureCompensator(makePtr<detail::BlocksGainCompensator>()); 571 stitcher->setBlender(makePtr<detail::MultiBandBlender>(try_use_gpu)); 572 573 return stitcher; 574 } 575 } // namespace cv 576