samples/cpp/stitching_detailed.cpp

   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 //
  42 //M*/
  43
  44 #include <iostream>
  45 #include <fstream>
  46 #include <string>
  47 #include "opencv2/opencv_modules.hpp"
  48 #include <opencv2/core/utility.hpp>
  49 #include "opencv2/imgcodecs.hpp"
  50 #include "opencv2/highgui.hpp"
  51 #include "opencv2/stitching/detail/autocalib.hpp"
  52 #include "opencv2/stitching/detail/blenders.hpp"
  53 #include "opencv2/stitching/detail/timelapsers.hpp"
  54 #include "opencv2/stitching/detail/camera.hpp"
  55 #include "opencv2/stitching/detail/exposure_compensate.hpp"
  56 #include "opencv2/stitching/detail/matchers.hpp"
  57 #include "opencv2/stitching/detail/motion_estimators.hpp"
  58 #include "opencv2/stitching/detail/seam_finders.hpp"
  59 #include "opencv2/stitching/detail/util.hpp"
  60 #include "opencv2/stitching/detail/warpers.hpp"
  61 #include "opencv2/stitching/warpers.hpp"
  62
  63 using namespace std;
  64 using namespace cv;
  65 using namespace cv::detail;
  66
  67 static void printUsage()
  68 {
  69     cout <<
  70         "Rotation model images stitcher.\n\n"
  71         "stitching_detailed img1 img2 [...imgN] [flags]\n\n"
  72         "Flags:\n"
  73         "  --preview\n"
  74         "      Run stitching in the preview mode. Works faster than usual mode,\n"
  75         "      but output image will have lower resolution.\n"
  76         "  --try_cuda (yes|no)\n"
  77         "      Try to use CUDA. The default value is 'no'. All default values\n"
  78         "      are for CPU mode.\n"
  79         "\nMotion Estimation Flags:\n"
  80         "  --work_megapix <float>\n"
  81         "      Resolution for image registration step. The default is 0.6 Mpx.\n"
  82         "  --features (surf|orb)\n"
  83         "      Type of features used for images matching. The default is surf.\n"
  84         "  --match_conf <float>\n"
  85         "      Confidence for feature matching step. The default is 0.65 for surf and 0.3 for orb.\n"
  86         "  --conf_thresh <float>\n"
  87         "      Threshold for two images are from the same panorama confidence.\n"
  88         "      The default is 1.0.\n"
  89         "  --ba (reproj|ray)\n"
  90         "      Bundle adjustment cost function. The default is ray.\n"
  91         "  --ba_refine_mask (mask)\n"
  92         "      Set refinement mask for bundle adjustment. It looks like 'x_xxx',\n"
  93         "      where 'x' means refine respective parameter and '_' means don't\n"
  94         "      refine one, and has the following format:\n"
  95         "      <fx><skew><ppx><aspect><ppy>. The default mask is 'xxxxx'. If bundle\n"
  96         "      adjustment doesn't support estimation of selected parameter then\n"
  97         "      the respective flag is ignored.\n"
  98         "  --wave_correct (no|horiz|vert)\n"
  99         "      Perform wave effect correction. The default is 'horiz'.\n"
 100         "  --save_graph <file_name>\n"
 101         "      Save matches graph represented in DOT language to <file_name> file.\n"
 102         "      Labels description: Nm is number of matches, Ni is number of inliers,\n"
 103         "      C is confidence.\n"
 104         "\nCompositing Flags:\n"
 105         "  --warp (plane|cylindrical|spherical|fisheye|stereographic|compressedPlaneA2B1|compressedPlaneA1.5B1|compressedPlanePortraitA2B1|compressedPlanePortraitA1.5B1|paniniA2B1|paniniA1.5B1|paniniPortraitA2B1|paniniPortraitA1.5B1|mercator|transverseMercator)\n"
 106         "      Warp surface type. The default is 'spherical'.\n"
 107         "  --seam_megapix <float>\n"
 108         "      Resolution for seam estimation step. The default is 0.1 Mpx.\n"
 109         "  --seam (no|voronoi|gc_color|gc_colorgrad)\n"
 110         "      Seam estimation method. The default is 'gc_color'.\n"
 111         "  --compose_megapix <float>\n"
 112         "      Resolution for compositing step. Use -1 for original resolution.\n"
 113         "      The default is -1.\n"
 114         "  --expos_comp (no|gain|gain_blocks)\n"
 115         "      Exposure compensation method. The default is 'gain_blocks'.\n"
 116         "  --blend (no|feather|multiband)\n"
 117         "      Blending method. The default is 'multiband'.\n"
 118         "  --blend_strength <float>\n"
 119         "      Blending strength from [0,100] range. The default is 5.\n"
 120         "  --output <result_img>\n"
 121         "      The default is 'result.jpg'.\n"
 122         "  --timelapse (as_is|crop) (range_width)\n"
 123         "      Output warped images separately as frames of a time lapse movie, with 'fixed_' prepended to input file names.\n";
 124 }
 125
 126
 127 // Default command line args
 128 vector<String> img_names;
 129 bool preview = false;
 130 bool try_cuda = false;
 131 double work_megapix = 0.6;
 132 double seam_megapix = 0.1;
 133 double compose_megapix = -1;
 134 float conf_thresh = 1.f;
 135 string features_type = "surf";
 136 string ba_cost_func = "ray";
 137 string ba_refine_mask = "xxxxx";
 138 bool do_wave_correct = true;
 139 WaveCorrectKind wave_correct = detail::WAVE_CORRECT_HORIZ;
 140 bool save_graph = false;
 141 std::string save_graph_to;
 142 string warp_type = "spherical";
 143 int expos_comp_type = ExposureCompensator::GAIN_BLOCKS;
 144 float match_conf = 0.3f;
 145 string seam_find_type = "gc_color";
 146 int blend_type = Blender::MULTI_BAND;
 147 int timelapse_type = Timelapser::AS_IS;
 148 float blend_strength = 5;
 149 string result_name = "result.jpg";
 150 bool timelapse = false;
 151 int timelapse_range = 5;
 152
 153
 154 static int parseCmdArgs(int argc, char** argv)
 155 {
 156     if (argc == 1)
 157     {
 158         printUsage();
 159         return -1;
 160     }
 161     for (int i = 1; i < argc; ++i)
 162     {
 163         if (string(argv[i]) == "--help" || string(argv[i]) == "/?")
 164         {
 165             printUsage();
 166             return -1;
 167         }
 168         else if (string(argv[i]) == "--preview")
 169         {
 170             preview = true;
 171         }
 172         else if (string(argv[i]) == "--try_cuda")
 173         {
 174             if (string(argv[i + 1]) == "no")
 175                 try_cuda = false;
 176             else if (string(argv[i + 1]) == "yes")
 177                 try_cuda = true;
 178             else
 179             {
 180                 cout << "Bad --try_cuda flag value\n";
 181                 return -1;
 182             }
 183             i++;
 184         }
 185         else if (string(argv[i]) == "--work_megapix")
 186         {
 187             work_megapix = atof(argv[i + 1]);
 188             i++;
 189         }
 190         else if (string(argv[i]) == "--seam_megapix")
 191         {
 192             seam_megapix = atof(argv[i + 1]);
 193             i++;
 194         }
 195         else if (string(argv[i]) == "--compose_megapix")
 196         {
 197             compose_megapix = atof(argv[i + 1]);
 198             i++;
 199         }
 200         else if (string(argv[i]) == "--result")
 201         {
 202             result_name = argv[i + 1];
 203             i++;
 204         }
 205         else if (string(argv[i]) == "--features")
 206         {
 207             features_type = argv[i + 1];
 208             if (features_type == "orb")
 209                 match_conf = 0.3f;
 210             i++;
 211         }
 212         else if (string(argv[i]) == "--match_conf")
 213         {
 214             match_conf = static_cast<float>(atof(argv[i + 1]));
 215             i++;
 216         }
 217         else if (string(argv[i]) == "--conf_thresh")
 218         {
 219             conf_thresh = static_cast<float>(atof(argv[i + 1]));
 220             i++;
 221         }
 222         else if (string(argv[i]) == "--ba")
 223         {
 224             ba_cost_func = argv[i + 1];
 225             i++;
 226         }
 227         else if (string(argv[i]) == "--ba_refine_mask")
 228         {
 229             ba_refine_mask = argv[i + 1];
 230             if (ba_refine_mask.size() != 5)
 231             {
 232                 cout << "Incorrect refinement mask length.\n";
 233                 return -1;
 234             }
 235             i++;
 236         }
 237         else if (string(argv[i]) == "--wave_correct")
 238         {
 239             if (string(argv[i + 1]) == "no")
 240                 do_wave_correct = false;
 241             else if (string(argv[i + 1]) == "horiz")
 242             {
 243                 do_wave_correct = true;
 244                 wave_correct = detail::WAVE_CORRECT_HORIZ;
 245             }
 246             else if (string(argv[i + 1]) == "vert")
 247             {
 248                 do_wave_correct = true;
 249                 wave_correct = detail::WAVE_CORRECT_VERT;
 250             }
 251             else
 252             {
 253                 cout << "Bad --wave_correct flag value\n";
 254                 return -1;
 255             }
 256             i++;
 257         }
 258         else if (string(argv[i]) == "--save_graph")
 259         {
 260             save_graph = true;
 261             save_graph_to = argv[i + 1];
 262             i++;
 263         }
 264         else if (string(argv[i]) == "--warp")
 265         {
 266             warp_type = string(argv[i + 1]);
 267             i++;
 268         }
 269         else if (string(argv[i]) == "--expos_comp")
 270         {
 271             if (string(argv[i + 1]) == "no")
 272                 expos_comp_type = ExposureCompensator::NO;
 273             else if (string(argv[i + 1]) == "gain")
 274                 expos_comp_type = ExposureCompensator::GAIN;
 275             else if (string(argv[i + 1]) == "gain_blocks")
 276                 expos_comp_type = ExposureCompensator::GAIN_BLOCKS;
 277             else
 278             {
 279                 cout << "Bad exposure compensation method\n";
 280                 return -1;
 281             }
 282             i++;
 283         }
 284         else if (string(argv[i]) == "--seam")
 285         {
 286             if (string(argv[i + 1]) == "no" ||
 287                 string(argv[i + 1]) == "voronoi" ||
 288                 string(argv[i + 1]) == "gc_color" ||
 289                 string(argv[i + 1]) == "gc_colorgrad" ||
 290                 string(argv[i + 1]) == "dp_color" ||
 291                 string(argv[i + 1]) == "dp_colorgrad")
 292                 seam_find_type = argv[i + 1];
 293             else
 294             {
 295                 cout << "Bad seam finding method\n";
 296                 return -1;
 297             }
 298             i++;
 299         }
 300         else if (string(argv[i]) == "--blend")
 301         {
 302             if (string(argv[i + 1]) == "no")
 303                 blend_type = Blender::NO;
 304             else if (string(argv[i + 1]) == "feather")
 305                 blend_type = Blender::FEATHER;
 306             else if (string(argv[i + 1]) == "multiband")
 307                 blend_type = Blender::MULTI_BAND;
 308             else
 309             {
 310                 cout << "Bad blending method\n";
 311                 return -1;
 312             }
 313             i++;
 314         }
 315         else if (string(argv[i]) == "--timelapse")
 316         {
 317             timelapse = true;
 318
 319             if (string(argv[i + 1]) == "as_is")
 320                 timelapse_type = Timelapser::AS_IS;
 321             else if (string(argv[i + 1]) == "crop")
 322                 timelapse_type = Timelapser::CROP;
 323             else
 324             {
 325                 cout << "Bad timelapse method\n";
 326                 return -1;
 327             }
 328             i++;
 329
 330             timelapse_range = atoi(argv[i + 1]);
 331             i++;
 332         }
 333         else if (string(argv[i]) == "--blend_strength")
 334         {
 335             blend_strength = static_cast<float>(atof(argv[i + 1]));
 336             i++;
 337         }
 338         else if (string(argv[i]) == "--output")
 339         {
 340             result_name = argv[i + 1];
 341             i++;
 342         }
 343         else
 344             img_names.push_back(argv[i]);
 345     }
 346     if (preview)
 347     {
 348         compose_megapix = 0.6;
 349     }
 350     return 0;
 351 }
 352
 353
 354 int main(int argc, char* argv[])
 355 {
 356 #if ENABLE_LOG
 357     int64 app_start_time = getTickCount();
 358 #endif
 359
 360 #if 0
 361     cv::setBreakOnError(true);
 362 #endif
 363
 364     int retval = parseCmdArgs(argc, argv);
 365     if (retval)
 366         return retval;
 367
 368     // Check if have enough images
 369     int num_images = static_cast<int>(img_names.size());
 370     if (num_images < 2)
 371     {
 372         LOGLN("Need more images");
 373         return -1;
 374     }
 375
 376     double work_scale = 1, seam_scale = 1, compose_scale = 1;
 377     bool is_work_scale_set = false, is_seam_scale_set = false, is_compose_scale_set = false;
 378
 379     LOGLN("Finding features...");
 380 #if ENABLE_LOG
 381     int64 t = getTickCount();
 382 #endif
 383
 384     Ptr<FeaturesFinder> finder;
 385     if (features_type == "surf")
 386     {
 387 #ifdef HAVE_OPENCV_NONFREE
 388         if (try_cuda && cuda::getCudaEnabledDeviceCount() > 0)
 389             finder = makePtr<SurfFeaturesFinderGpu>();
 390         else
 391 #endif
 392             finder = makePtr<SurfFeaturesFinder>();
 393     }
 394     else if (features_type == "orb")
 395     {
 396         finder = makePtr<OrbFeaturesFinder>();
 397     }
 398     else
 399     {
 400         cout << "Unknown 2D features type: '" << features_type << "'.\n";
 401         return -1;
 402     }
 403
 404     Mat full_img, img;
 405     vector<ImageFeatures> features(num_images);
 406     vector<Mat> images(num_images);
 407     vector<Size> full_img_sizes(num_images);
 408     double seam_work_aspect = 1;
 409
 410     for (int i = 0; i < num_images; ++i)
 411     {
 412         full_img = imread(img_names[i]);
 413         full_img_sizes[i] = full_img.size();
 414
 415         if (full_img.empty())
 416         {
 417             LOGLN("Can't open image " << img_names[i]);
 418             return -1;
 419         }
 420         if (work_megapix < 0)
 421         {
 422             img = full_img;
 423             work_scale = 1;
 424             is_work_scale_set = true;
 425         }
 426         else
 427         {
 428             if (!is_work_scale_set)
 429             {
 430                 work_scale = min(1.0, sqrt(work_megapix * 1e6 / full_img.size().area()));
 431                 is_work_scale_set = true;
 432             }
 433             resize(full_img, img, Size(), work_scale, work_scale);
 434         }
 435         if (!is_seam_scale_set)
 436         {
 437             seam_scale = min(1.0, sqrt(seam_megapix * 1e6 / full_img.size().area()));
 438             seam_work_aspect = seam_scale / work_scale;
 439             is_seam_scale_set = true;
 440         }
 441
 442         (*finder)(img, features[i]);
 443         features[i].img_idx = i;
 444         LOGLN("Features in image #" << i+1 << ": " << features[i].keypoints.size());
 445
 446         resize(full_img, img, Size(), seam_scale, seam_scale);
 447         images[i] = img.clone();
 448     }
 449
 450     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     vector<MatchesInfo> pairwise_matches;
 461     if (!timelapse)
 462     {
 463         BestOf2NearestMatcher matcher(try_cuda, match_conf);
 464         matcher(features, pairwise_matches);
 465         matcher.collectGarbage();
 466     }
 467     else
 468     {
 469         BestOf2NearestRangeMatcher matcher(timelapse_range, try_cuda, match_conf);
 470         matcher(features, pairwise_matches);
 471         matcher.collectGarbage();
 472     }
 473
 474     LOGLN("Pairwise matching, time: " << ((getTickCount() - t) / getTickFrequency()) << " sec");
 475
 476     // Check if we should save matches graph
 477     if (save_graph)
 478     {
 479         LOGLN("Saving matches graph...");
 480         ofstream f(save_graph_to.c_str());
 481         f << matchesGraphAsString(img_names, pairwise_matches, conf_thresh);
 482     }
 483
 484     // Leave only images we are sure are from the same panorama
 485     vector<int> indices = leaveBiggestComponent(features, pairwise_matches, conf_thresh);
 486     vector<Mat> img_subset;
 487     vector<String> img_names_subset;
 488     vector<Size> full_img_sizes_subset;
 489     for (size_t i = 0; i < indices.size(); ++i)
 490     {
 491         img_names_subset.push_back(img_names[indices[i]]);
 492         img_subset.push_back(images[indices[i]]);
 493         full_img_sizes_subset.push_back(full_img_sizes[indices[i]]);
 494     }
 495
 496     images = img_subset;
 497     img_names = img_names_subset;
 498     full_img_sizes = full_img_sizes_subset;
 499
 500     // Check if we still have enough images
 501     num_images = static_cast<int>(img_names.size());
 502     if (num_images < 2)
 503     {
 504         LOGLN("Need more images");
 505         return -1;
 506     }
 507
 508     HomographyBasedEstimator estimator;
 509     vector<CameraParams> cameras;
 510     if (!estimator(features, pairwise_matches, cameras))
 511     {
 512         cout << "Homography estimation failed.\n";
 513         return -1;
 514     }
 515
 516     for (size_t i = 0; i < cameras.size(); ++i)
 517     {
 518         Mat R;
 519         cameras[i].R.convertTo(R, CV_32F);
 520         cameras[i].R = R;
 521         LOGLN("Initial intrinsics #" << indices[i]+1 << ":\n" << cameras[i].K());
 522     }
 523
 524     Ptr<detail::BundleAdjusterBase> adjuster;
 525     if (ba_cost_func == "reproj") adjuster = makePtr<detail::BundleAdjusterReproj>();
 526     else if (ba_cost_func == "ray") adjuster = makePtr<detail::BundleAdjusterRay>();
 527     else
 528     {
 529         cout << "Unknown bundle adjustment cost function: '" << ba_cost_func << "'.\n";
 530         return -1;
 531     }
 532     adjuster->setConfThresh(conf_thresh);
 533     Mat_<uchar> refine_mask = Mat::zeros(3, 3, CV_8U);
 534     if (ba_refine_mask[0] == 'x') refine_mask(0,0) = 1;
 535     if (ba_refine_mask[1] == 'x') refine_mask(0,1) = 1;
 536     if (ba_refine_mask[2] == 'x') refine_mask(0,2) = 1;
 537     if (ba_refine_mask[3] == 'x') refine_mask(1,1) = 1;
 538     if (ba_refine_mask[4] == 'x') refine_mask(1,2) = 1;
 539     adjuster->setRefinementMask(refine_mask);
 540     if (!(*adjuster)(features, pairwise_matches, cameras))
 541     {
 542         cout << "Camera parameters adjusting failed.\n";
 543         return -1;
 544     }
 545
 546     // Find median focal length
 547
 548     vector<double> focals;
 549     for (size_t i = 0; i < cameras.size(); ++i)
 550     {
 551         LOGLN("Camera #" << indices[i]+1 << ":\n" << cameras[i].K());
 552         focals.push_back(cameras[i].focal);
 553     }
 554
 555     sort(focals.begin(), focals.end());
 556     float warped_image_scale;
 557     if (focals.size() % 2 == 1)
 558         warped_image_scale = static_cast<float>(focals[focals.size() / 2]);
 559     else
 560         warped_image_scale = static_cast<float>(focals[focals.size() / 2 - 1] + focals[focals.size() / 2]) * 0.5f;
 561
 562     if (do_wave_correct)
 563     {
 564         vector<Mat> rmats;
 565         for (size_t i = 0; i < cameras.size(); ++i)
 566             rmats.push_back(cameras[i].R.clone());
 567         waveCorrect(rmats, wave_correct);
 568         for (size_t i = 0; i < cameras.size(); ++i)
 569             cameras[i].R = rmats[i];
 570     }
 571
 572     LOGLN("Warping images (auxiliary)... ");
 573 #if ENABLE_LOG
 574     t = getTickCount();
 575 #endif
 576
 577     vector<Point> corners(num_images);
 578     vector<UMat> masks_warped(num_images);
 579     vector<UMat> images_warped(num_images);
 580     vector<Size> sizes(num_images);
 581     vector<UMat> masks(num_images);
 582
 583     // Preapre images masks
 584     for (int i = 0; i < num_images; ++i)
 585     {
 586         masks[i].create(images[i].size(), CV_8U);
 587         masks[i].setTo(Scalar::all(255));
 588     }
 589
 590     // Warp images and their masks
 591
 592     Ptr<WarperCreator> warper_creator;
 593 #ifdef HAVE_OPENCV_CUDAWARPING
 594     if (try_cuda && cuda::getCudaEnabledDeviceCount() > 0)
 595     {
 596         if (warp_type == "plane")
 597             warper_creator = makePtr<cv::PlaneWarperGpu>();
 598         else if (warp_type == "cylindrical")
 599             warper_creator = makePtr<cv::CylindricalWarperGpu>();
 600         else if (warp_type == "spherical")
 601             warper_creator = makePtr<cv::SphericalWarperGpu>();
 602     }
 603     else
 604 #endif
 605     {
 606         if (warp_type == "plane")
 607             warper_creator = makePtr<cv::PlaneWarper>();
 608         else if (warp_type == "cylindrical")
 609             warper_creator = makePtr<cv::CylindricalWarper>();
 610         else if (warp_type == "spherical")
 611             warper_creator = makePtr<cv::SphericalWarper>();
 612         else if (warp_type == "fisheye")
 613             warper_creator = makePtr<cv::FisheyeWarper>();
 614         else if (warp_type == "stereographic")
 615             warper_creator = makePtr<cv::StereographicWarper>();
 616         else if (warp_type == "compressedPlaneA2B1")
 617             warper_creator = makePtr<cv::CompressedRectilinearWarper>(2.0f, 1.0f);
 618         else if (warp_type == "compressedPlaneA1.5B1")
 619             warper_creator = makePtr<cv::CompressedRectilinearWarper>(1.5f, 1.0f);
 620         else if (warp_type == "compressedPlanePortraitA2B1")
 621             warper_creator = makePtr<cv::CompressedRectilinearPortraitWarper>(2.0f, 1.0f);
 622         else if (warp_type == "compressedPlanePortraitA1.5B1")
 623             warper_creator = makePtr<cv::CompressedRectilinearPortraitWarper>(1.5f, 1.0f);
 624         else if (warp_type == "paniniA2B1")
 625             warper_creator = makePtr<cv::PaniniWarper>(2.0f, 1.0f);
 626         else if (warp_type == "paniniA1.5B1")
 627             warper_creator = makePtr<cv::PaniniWarper>(1.5f, 1.0f);
 628         else if (warp_type == "paniniPortraitA2B1")
 629             warper_creator = makePtr<cv::PaniniPortraitWarper>(2.0f, 1.0f);
 630         else if (warp_type == "paniniPortraitA1.5B1")
 631             warper_creator = makePtr<cv::PaniniPortraitWarper>(1.5f, 1.0f);
 632         else if (warp_type == "mercator")
 633             warper_creator = makePtr<cv::MercatorWarper>();
 634         else if (warp_type == "transverseMercator")
 635             warper_creator = makePtr<cv::TransverseMercatorWarper>();
 636     }
 637
 638     if (!warper_creator)
 639     {
 640         cout << "Can't create the following warper '" << warp_type << "'\n";
 641         return 1;
 642     }
 643
 644     Ptr<RotationWarper> warper = warper_creator->create(static_cast<float>(warped_image_scale * seam_work_aspect));
 645
 646     for (int i = 0; i < num_images; ++i)
 647     {
 648         Mat_<float> K;
 649         cameras[i].K().convertTo(K, CV_32F);
 650         float swa = (float)seam_work_aspect;
 651         K(0,0) *= swa; K(0,2) *= swa;
 652         K(1,1) *= swa; K(1,2) *= swa;
 653
 654         corners[i] = warper->warp(images[i], K, cameras[i].R, INTER_LINEAR, BORDER_REFLECT, images_warped[i]);
 655         sizes[i] = images_warped[i].size();
 656
 657         warper->warp(masks[i], K, cameras[i].R, INTER_NEAREST, BORDER_CONSTANT, masks_warped[i]);
 658     }
 659
 660     vector<UMat> images_warped_f(num_images);
 661     for (int i = 0; i < num_images; ++i)
 662         images_warped[i].convertTo(images_warped_f[i], CV_32F);
 663
 664     LOGLN("Warping images, time: " << ((getTickCount() - t) / getTickFrequency()) << " sec");
 665
 666     Ptr<ExposureCompensator> compensator = ExposureCompensator::createDefault(expos_comp_type);
 667     compensator->feed(corners, images_warped, masks_warped);
 668
 669     Ptr<SeamFinder> seam_finder;
 670     if (seam_find_type == "no")
 671         seam_finder = makePtr<detail::NoSeamFinder>();
 672     else if (seam_find_type == "voronoi")
 673         seam_finder = makePtr<detail::VoronoiSeamFinder>();
 674     else if (seam_find_type == "gc_color")
 675     {
 676 #ifdef HAVE_OPENCV_CUDA
 677         if (try_cuda && cuda::getCudaEnabledDeviceCount() > 0)
 678             seam_finder = makePtr<detail::GraphCutSeamFinderGpu>(GraphCutSeamFinderBase::COST_COLOR);
 679         else
 680 #endif
 681             seam_finder = makePtr<detail::GraphCutSeamFinder>(GraphCutSeamFinderBase::COST_COLOR);
 682     }
 683     else if (seam_find_type == "gc_colorgrad")
 684     {
 685 #ifdef HAVE_OPENCV_CUDA
 686         if (try_cuda && cuda::getCudaEnabledDeviceCount() > 0)
 687             seam_finder = makePtr<detail::GraphCutSeamFinderGpu>(GraphCutSeamFinderBase::COST_COLOR_GRAD);
 688         else
 689 #endif
 690             seam_finder = makePtr<detail::GraphCutSeamFinder>(GraphCutSeamFinderBase::COST_COLOR_GRAD);
 691     }
 692     else if (seam_find_type == "dp_color")
 693         seam_finder = makePtr<detail::DpSeamFinder>(DpSeamFinder::COLOR);
 694     else if (seam_find_type == "dp_colorgrad")
 695         seam_finder = makePtr<detail::DpSeamFinder>(DpSeamFinder::COLOR_GRAD);
 696     if (!seam_finder)
 697     {
 698         cout << "Can't create the following seam finder '" << seam_find_type << "'\n";
 699         return 1;
 700     }
 701
 702     seam_finder->find(images_warped_f, corners, masks_warped);
 703
 704     // Release unused memory
 705     images.clear();
 706     images_warped.clear();
 707     images_warped_f.clear();
 708     masks.clear();
 709
 710     LOGLN("Compositing...");
 711 #if ENABLE_LOG
 712     t = getTickCount();
 713 #endif
 714
 715     Mat img_warped, img_warped_s;
 716     Mat dilated_mask, seam_mask, mask, mask_warped;
 717     Ptr<Blender> blender;
 718     Ptr<Timelapser> timelapser;
 719     //double compose_seam_aspect = 1;
 720     double compose_work_aspect = 1;
 721
 722     for (int img_idx = 0; img_idx < num_images; ++img_idx)
 723     {
 724         LOGLN("Compositing image #" << indices[img_idx]+1);
 725
 726         // Read image and resize it if necessary
 727         full_img = imread(img_names[img_idx]);
 728         if (!is_compose_scale_set)
 729         {
 730             if (compose_megapix > 0)
 731                 compose_scale = min(1.0, sqrt(compose_megapix * 1e6 / full_img.size().area()));
 732             is_compose_scale_set = true;
 733
 734             // Compute relative scales
 735             //compose_seam_aspect = compose_scale / seam_scale;
 736             compose_work_aspect = compose_scale / work_scale;
 737
 738             // Update warped image scale
 739             warped_image_scale *= static_cast<float>(compose_work_aspect);
 740             warper = warper_creator->create(warped_image_scale);
 741
 742             // Update corners and sizes
 743             for (int i = 0; i < num_images; ++i)
 744             {
 745                 // Update intrinsics
 746                 cameras[i].focal *= compose_work_aspect;
 747                 cameras[i].ppx *= compose_work_aspect;
 748                 cameras[i].ppy *= compose_work_aspect;
 749
 750                 // Update corner and size
 751                 Size sz = full_img_sizes[i];
 752                 if (std::abs(compose_scale - 1) > 1e-1)
 753                 {
 754                     sz.width = cvRound(full_img_sizes[i].width * compose_scale);
 755                     sz.height = cvRound(full_img_sizes[i].height * compose_scale);
 756                 }
 757
 758                 Mat K;
 759                 cameras[i].K().convertTo(K, CV_32F);
 760                 Rect roi = warper->warpRoi(sz, K, cameras[i].R);
 761                 corners[i] = roi.tl();
 762                 sizes[i] = roi.size();
 763             }
 764         }
 765         if (abs(compose_scale - 1) > 1e-1)
 766             resize(full_img, img, Size(), compose_scale, compose_scale);
 767         else
 768             img = full_img;
 769         full_img.release();
 770         Size img_size = img.size();
 771
 772         Mat K;
 773         cameras[img_idx].K().convertTo(K, CV_32F);
 774
 775         // Warp the current image
 776         warper->warp(img, K, cameras[img_idx].R, INTER_LINEAR, BORDER_REFLECT, img_warped);
 777
 778         // Warp the current image mask
 779         mask.create(img_size, CV_8U);
 780         mask.setTo(Scalar::all(255));
 781         warper->warp(mask, K, cameras[img_idx].R, INTER_NEAREST, BORDER_CONSTANT, mask_warped);
 782
 783         // Compensate exposure
 784         compensator->apply(img_idx, corners[img_idx], img_warped, mask_warped);
 785
 786         img_warped.convertTo(img_warped_s, CV_16S);
 787         img_warped.release();
 788         img.release();
 789         mask.release();
 790
 791         dilate(masks_warped[img_idx], dilated_mask, Mat());
 792         resize(dilated_mask, seam_mask, mask_warped.size());
 793         mask_warped = seam_mask & mask_warped;
 794
 795         if (!blender && !timelapse)
 796         {
 797             blender = Blender::createDefault(blend_type, try_cuda);
 798             Size dst_sz = resultRoi(corners, sizes).size();
 799             float blend_width = sqrt(static_cast<float>(dst_sz.area())) * blend_strength / 100.f;
 800             if (blend_width < 1.f)
 801                 blender = Blender::createDefault(Blender::NO, try_cuda);
 802             else if (blend_type == Blender::MULTI_BAND)
 803             {
 804                 MultiBandBlender* mb = dynamic_cast<MultiBandBlender*>(blender.get());
 805                 mb->setNumBands(static_cast<int>(ceil(log(blend_width)/log(2.)) - 1.));
 806                 LOGLN("Multi-band blender, number of bands: " << mb->numBands());
 807             }
 808             else if (blend_type == Blender::FEATHER)
 809             {
 810                 FeatherBlender* fb = dynamic_cast<FeatherBlender*>(blender.get());
 811                 fb->setSharpness(1.f/blend_width);
 812                 LOGLN("Feather blender, sharpness: " << fb->sharpness());
 813             }
 814             blender->prepare(corners, sizes);
 815         }
 816         else if (!timelapser)
 817         {
 818             CV_Assert(timelapse);
 819             timelapser = Timelapser::createDefault(timelapse_type);
 820             timelapser->initialize(corners, sizes);
 821         }
 822
 823         // Blend the current image
 824         if (timelapse)
 825         {
 826             timelapser->process(img_warped_s, Mat::ones(img_warped_s.size(), CV_8UC1), corners[img_idx]);
 827
 828             imwrite("fixed_" + img_names[img_idx], timelapser->getDst());
 829         }
 830         else
 831         {
 832             blender->feed(img_warped_s, mask_warped, corners[img_idx]);
 833         }
 834     }
 835
 836     if (!timelapse)
 837     {
 838         Mat result, result_mask;
 839         blender->blend(result, result_mask);
 840
 841         LOGLN("Compositing, time: " << ((getTickCount() - t) / getTickFrequency()) << " sec");
 842
 843         imwrite(result_name, result);
 844     }
 845
 846     LOGLN("Finished, total time: " << ((getTickCount() - app_start_time) / getTickFrequency()) << " sec");
 847     return 0;
 848 }