samples/cpp/stitching_detailed.cpp#
A detailed example on image stitching
1#include <iostream>
2#include <fstream>
3#include <string>
4#include "opencv2/opencv_modules.hpp"
5#include <opencv2/core/utility.hpp>
6#include "opencv2/imgcodecs.hpp"
7#include "opencv2/highgui.hpp"
8#include "opencv2/stitching/detail/autocalib.hpp"
9#include "opencv2/stitching/detail/blenders.hpp"
10#include "opencv2/stitching/detail/timelapsers.hpp"
11#include "opencv2/stitching/detail/camera.hpp"
12#include "opencv2/stitching/detail/exposure_compensate.hpp"
13#include "opencv2/stitching/detail/matchers.hpp"
14#include "opencv2/stitching/detail/motion_estimators.hpp"
15#include "opencv2/stitching/detail/seam_finders.hpp"
16#include "opencv2/stitching/detail/warpers.hpp"
17#include "opencv2/stitching/warpers.hpp"
18#include "opencv2/features.hpp"
19#include "opencv2/core/ocl.hpp"
20
21#ifdef HAVE_OPENCV_XFEATURES2D
22#include "opencv2/xfeatures2d.hpp"
23#include "opencv2/xfeatures2d/nonfree.hpp"
24#endif
25
26#define ENABLE_LOG 1
27#define LOG(msg) std::cout << msg
28#define LOGLN(msg) std::cout << msg << std::endl
29
30using namespace std;
31using namespace cv;
32using namespace cv::detail;
33
34static void printUsage(char** argv)
35{
36 cout <<
37 "Rotation model images stitcher.\n\n"
38 << argv[0] << " img1 img2 [...imgN] [flags]\n\n"
39 "Flags:\n"
40 " --preview\n"
41 " Run stitching in the preview mode. Works faster than usual mode,\n"
42 " but output image will have lower resolution.\n"
43 " --try_cuda (yes|no)\n"
44 " Try to use CUDA. The default value is 'no'. All default values\n"
45 " are for CPU mode.\n"
46 "\nMotion Estimation Flags:\n"
47 " --work_megapix <float>\n"
48 " Resolution for image registration step. The default is 0.6 Mpx.\n"
49 " --features (surf|orb|sift|akaze|aliked)\n"
50 " Type of features used for images matching.\n"
51 " The default is surf if available, orb otherwise.\n"
52 " When using 'aliked', requires --matcher lightglue and DNN model paths.\n"
53 " --matcher (homography|affine)\n"
54 " Matcher used for pairwise image matching.\n"
55 " --estimator (homography|affine)\n"
56 " Type of estimator used for transformation estimation.\n"
57 " --match_conf <float>\n"
58 " Confidence for feature matching step. The default is 0.65 for surf and 0.3 for orb.\n"
59 " --conf_thresh <float>\n"
60 " Threshold for two images are from the same panorama confidence.\n"
61 " The default is 1.0.\n"
62 " --ba (no|reproj|ray|affine)\n"
63 " Bundle adjustment cost function. The default is ray.\n"
64 " --ba_refine_mask (mask)\n"
65 " Set refinement mask for bundle adjustment. It looks like 'x_xxx',\n"
66 " where 'x' means refine respective parameter and '_' means don't\n"
67 " refine one, and has the following format:\n"
68 " <fx><skew><ppx><aspect><ppy>. The default mask is 'xxxxx'. If bundle\n"
69 " adjustment doesn't support estimation of selected parameter then\n"
70 " the respective flag is ignored.\n"
71 " --wave_correct (no|horiz|vert)\n"
72 " Perform wave effect correction. The default is 'horiz'.\n"
73 " --save_graph <file_name>\n"
74 " Save matches graph represented in DOT language to <file_name> file.\n"
75 " Labels description: Nm is number of matches, Ni is number of inliers,\n"
76 " C is confidence.\n"
77 "\nCompositing Flags:\n"
78 " --warp (affine|plane|cylindrical|spherical|fisheye|stereographic|compressedPlaneA2B1|compressedPlaneA1.5B1|compressedPlanePortraitA2B1|compressedPlanePortraitA1.5B1|paniniA2B1|paniniA1.5B1|paniniPortraitA2B1|paniniPortraitA1.5B1|mercator|transverseMercator)\n"
79 " Warp surface type. The default is 'spherical'.\n"
80 " --seam_megapix <float>\n"
81 " Resolution for seam estimation step. The default is 0.1 Mpx.\n"
82 " --seam (no|voronoi|gc_color|gc_colorgrad)\n"
83 " Seam estimation method. The default is 'gc_color'.\n"
84 " --compose_megapix <float>\n"
85 " Resolution for compositing step. Use -1 for original resolution.\n"
86 " The default is -1.\n"
87 " --expos_comp (no|gain|gain_blocks|channels|channels_blocks)\n"
88 " Exposure compensation method. The default is 'gain_blocks'.\n"
89 " --expos_comp_nr_feeds <int>\n"
90 " Number of exposure compensation feed. The default is 1.\n"
91 " --expos_comp_nr_filtering <int>\n"
92 " Number of filtering iterations of the exposure compensation gains.\n"
93 " Only used when using a block exposure compensation method.\n"
94 " The default is 2.\n"
95 " --expos_comp_block_size <int>\n"
96 " BLock size in pixels used by the exposure compensator.\n"
97 " Only used when using a block exposure compensation method.\n"
98 " The default is 32.\n"
99 " --blend (no|feather|multiband)\n"
100 " Blending method. The default is 'multiband'.\n"
101 " --blend_strength <float>\n"
102 " Blending strength from [0,100] range. The default is 5.\n"
103 " --output <result_img>\n"
104 " The default is 'result.jpg'.\n"
105 " --timelapse (as_is|crop) \n"
106 " Output warped images separately as frames of a time lapse movie, with 'fixed_' prepended to input file names.\n"
107 " --rangewidth <int>\n"
108 " uses range_width to limit number of images to match with.\n"
109 "\nDNN Feature Options (when --features aliked --matcher lightglue):\n"
110 " --aliked_model <path>\n"
111 " Path to ALIKED ONNX model file.\n"
112 " --lightglue_model <path>\n"
113 " Path to LightGlue ONNX model file (for ALIKED descriptors).\n"
114 " --lg_score_thresh <float>\n"
115 " LightGlue confidence threshold. The default is 0.0 (accept all).\n";
116}
117
118
119// Default command line args
120vector<String> img_names;
121bool preview = false;
122bool try_cuda = false;
123double work_megapix = 0.6;
124double seam_megapix = 0.1;
125double compose_megapix = -1;
126float conf_thresh = 1.f;
127#ifdef HAVE_OPENCV_XFEATURES2D
128string features_type = "surf";
129float match_conf = 0.65f;
130#else
131string features_type = "orb";
132float match_conf = 0.3f;
133#endif
134string matcher_type = "homography";
135string estimator_type = "homography";
136string ba_cost_func = "ray";
137string ba_refine_mask = "xxxxx";
138bool do_wave_correct = true;
139WaveCorrectKind wave_correct = detail::WAVE_CORRECT_HORIZ;
140bool save_graph = false;
141std::string save_graph_to;
142string warp_type = "spherical";
143int expos_comp_type = ExposureCompensator::GAIN_BLOCKS;
144int expos_comp_nr_feeds = 1;
145int expos_comp_nr_filtering = 2;
146int expos_comp_block_size = 32;
147string seam_find_type = "gc_color";
148int blend_type = Blender::MULTI_BAND;
149int timelapse_type = Timelapser::AS_IS;
150float blend_strength = 5;
151string result_name = "result.jpg";
152bool timelapse = false;
153int range_width = -1;
154String aliked_model_path;
155String lightglue_model_path;
156float lg_score_thresh = 0.0f;
157
158
159static int parseCmdArgs(int argc, char** argv)
160{
161 if (argc == 1)
162 {
163 printUsage(argv);
164 return -1;
165 }
166 for (int i = 1; i < argc; ++i)
167 {
168 if (string(argv[i]) == "--help" || string(argv[i]) == "/?")
169 {
170 printUsage(argv);
171 return -1;
172 }
173 else if (string(argv[i]) == "--preview")
174 {
175 preview = true;
176 }
177 else if (string(argv[i]) == "--try_cuda")
178 {
179 if (string(argv[i + 1]) == "no")
180 try_cuda = false;
181 else if (string(argv[i + 1]) == "yes")
182 try_cuda = true;
183 else
184 {
185 cout << "Bad --try_cuda flag value\n";
186 return -1;
187 }
188 i++;
189 }
190 else if (string(argv[i]) == "--work_megapix")
191 {
192 work_megapix = atof(argv[i + 1]);
193 i++;
194 }
195 else if (string(argv[i]) == "--seam_megapix")
196 {
197 seam_megapix = atof(argv[i + 1]);
198 i++;
199 }
200 else if (string(argv[i]) == "--compose_megapix")
201 {
202 compose_megapix = atof(argv[i + 1]);
203 i++;
204 }
205 else if (string(argv[i]) == "--result")
206 {
207 result_name = argv[i + 1];
208 i++;
209 }
210 else if (string(argv[i]) == "--features")
211 {
212 features_type = argv[i + 1];
213 if (string(features_type) == "orb")
214 match_conf = 0.3f;
215 i++;
216 }
217 else if (string(argv[i]) == "--matcher")
218 {
219 if (string(argv[i + 1]) == "homography" || string(argv[i + 1]) == "affine" || string(argv[i + 1]) == "lightglue")
220 matcher_type = argv[i + 1];
221 else
222 {
223 cout << "Bad --matcher flag value\n";
224 return -1;
225 }
226 i++;
227 }
228 else if (string(argv[i]) == "--estimator")
229 {
230 if (string(argv[i + 1]) == "homography" || string(argv[i + 1]) == "affine")
231 estimator_type = argv[i + 1];
232 else
233 {
234 cout << "Bad --estimator flag value\n";
235 return -1;
236 }
237 i++;
238 }
239 else if (string(argv[i]) == "--match_conf")
240 {
241 match_conf = static_cast<float>(atof(argv[i + 1]));
242 i++;
243 }
244 else if (string(argv[i]) == "--conf_thresh")
245 {
246 conf_thresh = static_cast<float>(atof(argv[i + 1]));
247 i++;
248 }
249 else if (string(argv[i]) == "--ba")
250 {
251 ba_cost_func = argv[i + 1];
252 i++;
253 }
254 else if (string(argv[i]) == "--ba_refine_mask")
255 {
256 ba_refine_mask = argv[i + 1];
257 if (ba_refine_mask.size() != 5)
258 {
259 cout << "Incorrect refinement mask length.\n";
260 return -1;
261 }
262 i++;
263 }
264 else if (string(argv[i]) == "--wave_correct")
265 {
266 if (string(argv[i + 1]) == "no")
267 do_wave_correct = false;
268 else if (string(argv[i + 1]) == "horiz")
269 {
270 do_wave_correct = true;
271 wave_correct = detail::WAVE_CORRECT_HORIZ;
272 }
273 else if (string(argv[i + 1]) == "vert")
274 {
275 do_wave_correct = true;
276 wave_correct = detail::WAVE_CORRECT_VERT;
277 }
278 else
279 {
280 cout << "Bad --wave_correct flag value\n";
281 return -1;
282 }
283 i++;
284 }
285 else if (string(argv[i]) == "--save_graph")
286 {
287 save_graph = true;
288 save_graph_to = argv[i + 1];
289 i++;
290 }
291 else if (string(argv[i]) == "--warp")
292 {
293 warp_type = string(argv[i + 1]);
294 i++;
295 }
296 else if (string(argv[i]) == "--expos_comp")
297 {
298 if (string(argv[i + 1]) == "no")
299 expos_comp_type = ExposureCompensator::NO;
300 else if (string(argv[i + 1]) == "gain")
301 expos_comp_type = ExposureCompensator::GAIN;
302 else if (string(argv[i + 1]) == "gain_blocks")
303 expos_comp_type = ExposureCompensator::GAIN_BLOCKS;
304 else if (string(argv[i + 1]) == "channels")
305 expos_comp_type = ExposureCompensator::CHANNELS;
306 else if (string(argv[i + 1]) == "channels_blocks")
307 expos_comp_type = ExposureCompensator::CHANNELS_BLOCKS;
308 else
309 {
310 cout << "Bad exposure compensation method\n";
311 return -1;
312 }
313 i++;
314 }
315 else if (string(argv[i]) == "--expos_comp_nr_feeds")
316 {
317 expos_comp_nr_feeds = atoi(argv[i + 1]);
318 i++;
319 }
320 else if (string(argv[i]) == "--expos_comp_nr_filtering")
321 {
322 expos_comp_nr_filtering = atoi(argv[i + 1]);
323 i++;
324 }
325 else if (string(argv[i]) == "--expos_comp_block_size")
326 {
327 expos_comp_block_size = atoi(argv[i + 1]);
328 i++;
329 }
330 else if (string(argv[i]) == "--seam")
331 {
332 if (string(argv[i + 1]) == "no" ||
333 string(argv[i + 1]) == "voronoi" ||
334 string(argv[i + 1]) == "gc_color" ||
335 string(argv[i + 1]) == "gc_colorgrad" ||
336 string(argv[i + 1]) == "dp_color" ||
337 string(argv[i + 1]) == "dp_colorgrad")
338 seam_find_type = argv[i + 1];
339 else
340 {
341 cout << "Bad seam finding method\n";
342 return -1;
343 }
344 i++;
345 }
346 else if (string(argv[i]) == "--blend")
347 {
348 if (string(argv[i + 1]) == "no")
349 blend_type = Blender::NO;
350 else if (string(argv[i + 1]) == "feather")
351 blend_type = Blender::FEATHER;
352 else if (string(argv[i + 1]) == "multiband")
353 blend_type = Blender::MULTI_BAND;
354 else
355 {
356 cout << "Bad blending method\n";
357 return -1;
358 }
359 i++;
360 }
361 else if (string(argv[i]) == "--timelapse")
362 {
363 timelapse = true;
364
365 if (string(argv[i + 1]) == "as_is")
366 timelapse_type = Timelapser::AS_IS;
367 else if (string(argv[i + 1]) == "crop")
368 timelapse_type = Timelapser::CROP;
369 else
370 {
371 cout << "Bad timelapse method\n";
372 return -1;
373 }
374 i++;
375 }
376 else if (string(argv[i]) == "--rangewidth")
377 {
378 range_width = atoi(argv[i + 1]);
379 i++;
380 }
381 else if (string(argv[i]) == "--blend_strength")
382 {
383 blend_strength = static_cast<float>(atof(argv[i + 1]));
384 i++;
385 }
386 else if (string(argv[i]) == "--output")
387 {
388 result_name = argv[i + 1];
389 i++;
390 }
391 else if (string(argv[i]) == "--aliked_model")
392 {
393 aliked_model_path = argv[i + 1];
394 i++;
395 }
396 else if (string(argv[i]) == "--lightglue_model")
397 {
398 lightglue_model_path = argv[i + 1];
399 i++;
400 }
401 else if (string(argv[i]) == "--lg_score_thresh")
402 {
403 lg_score_thresh = static_cast<float>(atof(argv[i + 1]));
404 i++;
405 }
406 else
407 img_names.push_back(argv[i]);
408 }
409 if (preview)
410 {
411 compose_megapix = 0.6;
412 }
413
414 // Validate DNN options
415 if (features_type == "aliked" && matcher_type != "lightglue")
416 {
417 cout << "Error: --features aliked requires --matcher lightglue\n";
418 return -1;
419 }
420 if (features_type == "aliked" && (aliked_model_path.empty() || lightglue_model_path.empty()))
421 {
422 cout << "Error: --features aliked requires --aliked_model and --lightglue_model\n";
423 return -1;
424 }
425
426 return 0;
427}
428
429
430int main(int argc, char* argv[])
431{
432#if ENABLE_LOG
433 int64 app_start_time = getTickCount();
434#endif
435
436#if 0
437 cv::setBreakOnError(true);
438#endif
439
440 int retval = parseCmdArgs(argc, argv);
441 if (retval)
442 return retval;
443
444 // Disable OpenCL for DNN-based features to avoid backend sync issues
445 bool use_aliked = (features_type == "aliked");
446 if (use_aliked)
447 cv::ocl::setUseOpenCL(false);
448
449 // Check if have enough images
450 int num_images = static_cast<int>(img_names.size());
451 if (num_images < 2)
452 {
453 LOGLN("Need more images");
454 return -1;
455 }
456
457 double work_scale = 1, seam_scale = 1, compose_scale = 1;
458 bool is_work_scale_set = false, is_seam_scale_set = false, is_compose_scale_set = false;
459
460 LOGLN("Finding features...");
461#if ENABLE_LOG
462 int64 t = getTickCount();
463#endif
464
465 Ptr<Feature2D> finder;
466 if (use_aliked)
467 {
468 // ALIKED will be created per-image in the loop below
469 }
470 else if (features_type == "orb")
471 {
472 finder = ORB::create();
473 }
474 else if (features_type == "akaze")
475 {
476#ifdef HAVE_OPENCV_XFEATURES2D
477 finder = xfeatures2d::AKAZE::create();
478#else
479 cout << "OpenCV is built without opencv_contrib modules. AKAZE algorithm is not available!" << std::endl;
480 return -1;
481#endif
482 }
483 else if (features_type == "surf")
484 {
485#if defined(HAVE_OPENCV_XFEATURES2D) && defined(HAVE_OPENCV_NONFREE)
486 finder = xfeatures2d::SURF::create();
487#else
488 cout << "OpenCV is built without NONFREE modules. SURF algorithm is not available!" << std::endl;
489 return -1;
490#endif
491 }
492 else if (features_type == "sift")
493 {
494 finder = SIFT::create();
495 }
496 else
497 {
498 cout << "Unknown 2D features type: '" << features_type << "'.\n";
499 return -1;
500 }
501
502 Mat full_img, img;
503 vector<ImageFeatures> features(num_images);
504 vector<Mat> images(num_images);
505 vector<Size> full_img_sizes(num_images);
506 double seam_work_aspect = 1;
507
508 for (int i = 0; i < num_images; ++i)
509 {
510 full_img = imread(samples::findFile(img_names[i]));
511 full_img_sizes[i] = full_img.size();
512
513 if (full_img.empty())
514 {
515 LOGLN("Can't open image " << img_names[i]);
516 return -1;
517 }
518 if (work_megapix < 0)
519 {
520 img = full_img;
521 work_scale = 1;
522 is_work_scale_set = true;
523 }
524 else
525 {
526 if (!is_work_scale_set)
527 {
528 work_scale = min(1.0, sqrt(work_megapix * 1e6 / full_img.size().area()));
529 is_work_scale_set = true;
530 }
531 resize(full_img, img, Size(), work_scale, work_scale, INTER_LINEAR_EXACT);
532 }
533 if (!is_seam_scale_set)
534 {
535 seam_scale = min(1.0, sqrt(seam_megapix * 1e6 / full_img.size().area()));
536 seam_work_aspect = seam_scale / work_scale;
537 is_seam_scale_set = true;
538 }
539
540 if (use_aliked)
541 {
542 Ptr<ALIKED> aliked = ALIKED::create(aliked_model_path);
543 computeImageFeatures(aliked, img, features[i]);
544 }
545 else
546 {
547 computeImageFeatures(finder, img, features[i]);
548 }
549 features[i].img_idx = i;
550 LOGLN("Features in image #" << i+1 << ": " << features[i].keypoints.size());
551
552 resize(full_img, img, Size(), seam_scale, seam_scale, INTER_LINEAR_EXACT);
553 images[i] = img.clone();
554 }
555
556 full_img.release();
557 img.release();
558
559 LOGLN("Finding features, time: " << ((getTickCount() - t) / getTickFrequency()) << " sec");
560
561 LOG("Pairwise matching");
562#if ENABLE_LOG
563 t = getTickCount();
564#endif
565 vector<MatchesInfo> pairwise_matches;
566 Ptr<FeaturesMatcher> matcher;
567 if (use_aliked && matcher_type == "lightglue")
568 {
569 Ptr<LightGlueMatcher> lg = LightGlueMatcher::create(lightglue_model_path);
570 Ptr<LightGlueFeaturesMatcher> lgMatcher = makePtr<LightGlueFeaturesMatcher>(lg);
571 lgMatcher->setScoreThreshold(lg_score_thresh);
572 matcher = lgMatcher;
573 }
574 else if (matcher_type == "affine")
575 matcher = makePtr<AffineBestOf2NearestMatcher>(false, try_cuda, match_conf);
576 else if (range_width==-1)
577 matcher = makePtr<BestOf2NearestMatcher>(try_cuda, match_conf);
578 else
579 matcher = makePtr<BestOf2NearestRangeMatcher>(range_width, try_cuda, match_conf);
580
581 (*matcher)(features, pairwise_matches);
582 matcher->collectGarbage();
583
584 LOGLN("Pairwise matching, time: " << ((getTickCount() - t) / getTickFrequency()) << " sec");
585
586 // Check if we should save matches graph
587 if (save_graph)
588 {
589 LOGLN("Saving matches graph...");
590 ofstream f(save_graph_to.c_str());
591 f << matchesGraphAsString(img_names, pairwise_matches, conf_thresh);
592 }
593
594 // Leave only images we are sure are from the same panorama
595 vector<int> indices = leaveBiggestComponent(features, pairwise_matches, conf_thresh);
596 vector<Mat> img_subset;
597 vector<String> img_names_subset;
598 vector<Size> full_img_sizes_subset;
599 for (size_t i = 0; i < indices.size(); ++i)
600 {
601 img_names_subset.push_back(img_names[indices[i]]);
602 img_subset.push_back(images[indices[i]]);
603 full_img_sizes_subset.push_back(full_img_sizes[indices[i]]);
604 }
605
606 images = img_subset;
607 img_names = img_names_subset;
608 full_img_sizes = full_img_sizes_subset;
609
610 // Check if we still have enough images
611 num_images = static_cast<int>(img_names.size());
612 if (num_images < 2)
613 {
614 LOGLN("Need more images");
615 return -1;
616 }
617
618 Ptr<Estimator> estimator;
619 if (estimator_type == "affine")
620 estimator = makePtr<AffineBasedEstimator>();
621 else
622 estimator = makePtr<HomographyBasedEstimator>();
623
624 vector<CameraParams> cameras;
625 if (!(*estimator)(features, pairwise_matches, cameras))
626 {
627 cout << "Homography estimation failed.\n";
628 return -1;
629 }
630
631 for (size_t i = 0; i < cameras.size(); ++i)
632 {
633 Mat R;
634 cameras[i].R.convertTo(R, CV_32F);
635 cameras[i].R = R;
636 LOGLN("Initial camera intrinsics #" << indices[i]+1 << ":\nK:\n" << cameras[i].K() << "\nR:\n" << cameras[i].R);
637 }
638
639 Ptr<detail::BundleAdjusterBase> adjuster;
640 if (ba_cost_func == "reproj") adjuster = makePtr<detail::BundleAdjusterReproj>();
641 else if (ba_cost_func == "ray") adjuster = makePtr<detail::BundleAdjusterRay>();
642 else if (ba_cost_func == "affine") adjuster = makePtr<detail::BundleAdjusterAffinePartial>();
643 else if (ba_cost_func == "no") adjuster = makePtr<NoBundleAdjuster>();
644 else
645 {
646 cout << "Unknown bundle adjustment cost function: '" << ba_cost_func << "'.\n";
647 return -1;
648 }
649 adjuster->setConfThresh(conf_thresh);
650 Mat_<uchar> refine_mask = Mat::zeros(3, 3, CV_8U);
651 if (ba_refine_mask[0] == 'x') refine_mask(0,0) = 1;
652 if (ba_refine_mask[1] == 'x') refine_mask(0,1) = 1;
653 if (ba_refine_mask[2] == 'x') refine_mask(0,2) = 1;
654 if (ba_refine_mask[3] == 'x') refine_mask(1,1) = 1;
655 if (ba_refine_mask[4] == 'x') refine_mask(1,2) = 1;
656 adjuster->setRefinementMask(refine_mask);
657 if (!(*adjuster)(features, pairwise_matches, cameras))
658 {
659 cout << "Camera parameters adjusting failed.\n";
660 return -1;
661 }
662
663 // Find median focal length
664
665 vector<double> focals;
666 for (size_t i = 0; i < cameras.size(); ++i)
667 {
668 LOGLN("Camera #" << indices[i]+1 << ":\nK:\n" << cameras[i].K() << "\nR:\n" << cameras[i].R);
669 focals.push_back(cameras[i].focal);
670 }
671
672 sort(focals.begin(), focals.end());
673 float warped_image_scale;
674 if (focals.size() % 2 == 1)
675 warped_image_scale = static_cast<float>(focals[focals.size() / 2]);
676 else
677 warped_image_scale = static_cast<float>(focals[focals.size() / 2 - 1] + focals[focals.size() / 2]) * 0.5f;
678
679 if (do_wave_correct)
680 {
681 vector<Mat> rmats;
682 for (size_t i = 0; i < cameras.size(); ++i)
683 rmats.push_back(cameras[i].R.clone());
684 waveCorrect(rmats, wave_correct);
685 for (size_t i = 0; i < cameras.size(); ++i)
686 cameras[i].R = rmats[i];
687 }
688
689 LOGLN("Warping images (auxiliary)... ");
690#if ENABLE_LOG
691 t = getTickCount();
692#endif
693
694 vector<Point> corners(num_images);
695 vector<UMat> masks_warped(num_images);
696 vector<UMat> images_warped(num_images);
697 vector<Size> sizes(num_images);
698 vector<UMat> masks(num_images);
699
700 // Prepare images masks
701 for (int i = 0; i < num_images; ++i)
702 {
703 masks[i].create(images[i].size(), CV_8U);
704 masks[i].setTo(Scalar::all(255));
705 }
706
707 // Warp images and their masks
708
709 Ptr<WarperCreator> warper_creator;
710#ifdef HAVE_OPENCV_CUDAWARPING
711 if (try_cuda && cuda::getCudaEnabledDeviceCount() > 0)
712 {
713 if (warp_type == "plane")
714 warper_creator = makePtr<cv::PlaneWarperGpu>();
715 else if (warp_type == "cylindrical")
716 warper_creator = makePtr<cv::CylindricalWarperGpu>();
717 else if (warp_type == "spherical")
718 warper_creator = makePtr<cv::SphericalWarperGpu>();
719 }
720 else
721#endif
722 {
723 if (warp_type == "plane")
724 warper_creator = makePtr<cv::PlaneWarper>();
725 else if (warp_type == "affine")
726 warper_creator = makePtr<cv::AffineWarper>();
727 else if (warp_type == "cylindrical")
728 warper_creator = makePtr<cv::CylindricalWarper>();
729 else if (warp_type == "spherical")
730 warper_creator = makePtr<cv::SphericalWarper>();
731 else if (warp_type == "fisheye")
732 warper_creator = makePtr<cv::FisheyeWarper>();
733 else if (warp_type == "stereographic")
734 warper_creator = makePtr<cv::StereographicWarper>();
735 else if (warp_type == "compressedPlaneA2B1")
736 warper_creator = makePtr<cv::CompressedRectilinearWarper>(2.0f, 1.0f);
737 else if (warp_type == "compressedPlaneA1.5B1")
738 warper_creator = makePtr<cv::CompressedRectilinearWarper>(1.5f, 1.0f);
739 else if (warp_type == "compressedPlanePortraitA2B1")
740 warper_creator = makePtr<cv::CompressedRectilinearPortraitWarper>(2.0f, 1.0f);
741 else if (warp_type == "compressedPlanePortraitA1.5B1")
742 warper_creator = makePtr<cv::CompressedRectilinearPortraitWarper>(1.5f, 1.0f);
743 else if (warp_type == "paniniA2B1")
744 warper_creator = makePtr<cv::PaniniWarper>(2.0f, 1.0f);
745 else if (warp_type == "paniniA1.5B1")
746 warper_creator = makePtr<cv::PaniniWarper>(1.5f, 1.0f);
747 else if (warp_type == "paniniPortraitA2B1")
748 warper_creator = makePtr<cv::PaniniPortraitWarper>(2.0f, 1.0f);
749 else if (warp_type == "paniniPortraitA1.5B1")
750 warper_creator = makePtr<cv::PaniniPortraitWarper>(1.5f, 1.0f);
751 else if (warp_type == "mercator")
752 warper_creator = makePtr<cv::MercatorWarper>();
753 else if (warp_type == "transverseMercator")
754 warper_creator = makePtr<cv::TransverseMercatorWarper>();
755 }
756
757 if (!warper_creator)
758 {
759 cout << "Can't create the following warper '" << warp_type << "'\n";
760 return 1;
761 }
762
763 Ptr<RotationWarper> warper = warper_creator->create(static_cast<float>(warped_image_scale * seam_work_aspect));
764
765 for (int i = 0; i < num_images; ++i)
766 {
767 Mat_<float> K;
768 cameras[i].K().convertTo(K, CV_32F);
769 float swa = (float)seam_work_aspect;
770 K(0,0) *= swa; K(0,2) *= swa;
771 K(1,1) *= swa; K(1,2) *= swa;
772
773 corners[i] = warper->warp(images[i], K, cameras[i].R, INTER_LINEAR, BORDER_REFLECT, images_warped[i]);
774 sizes[i] = images_warped[i].size();
775
776 warper->warp(masks[i], K, cameras[i].R, INTER_NEAREST, BORDER_CONSTANT, masks_warped[i]);
777 }
778
779 vector<UMat> images_warped_f(num_images);
780 for (int i = 0; i < num_images; ++i)
781 images_warped[i].convertTo(images_warped_f[i], CV_32F);
782
783 LOGLN("Warping images, time: " << ((getTickCount() - t) / getTickFrequency()) << " sec");
784
785 LOGLN("Compensating exposure...");
786#if ENABLE_LOG
787 t = getTickCount();
788#endif
789
790 Ptr<ExposureCompensator> compensator = ExposureCompensator::createDefault(expos_comp_type);
791 if (dynamic_cast<GainCompensator*>(compensator.get()))
792 {
793 GainCompensator* gcompensator = dynamic_cast<GainCompensator*>(compensator.get());
794 gcompensator->setNrFeeds(expos_comp_nr_feeds);
795 }
796
797 if (dynamic_cast<ChannelsCompensator*>(compensator.get()))
798 {
799 ChannelsCompensator* ccompensator = dynamic_cast<ChannelsCompensator*>(compensator.get());
800 ccompensator->setNrFeeds(expos_comp_nr_feeds);
801 }
802
803 if (dynamic_cast<BlocksCompensator*>(compensator.get()))
804 {
805 BlocksCompensator* bcompensator = dynamic_cast<BlocksCompensator*>(compensator.get());
806 bcompensator->setNrFeeds(expos_comp_nr_feeds);
807 bcompensator->setNrGainsFilteringIterations(expos_comp_nr_filtering);
808 bcompensator->setBlockSize(expos_comp_block_size, expos_comp_block_size);
809 }
810
811 compensator->feed(corners, images_warped, masks_warped);
812
813 LOGLN("Compensating exposure, time: " << ((getTickCount() - t) / getTickFrequency()) << " sec");
814
815 LOGLN("Finding seams...");
816#if ENABLE_LOG
817 t = getTickCount();
818#endif
819
820 Ptr<SeamFinder> seam_finder;
821 if (seam_find_type == "no")
822 seam_finder = makePtr<detail::NoSeamFinder>();
823 else if (seam_find_type == "voronoi")
824 seam_finder = makePtr<detail::VoronoiSeamFinder>();
825 else if (seam_find_type == "gc_color")
826 {
827#ifdef HAVE_OPENCV_CUDALEGACY
828 if (try_cuda && cuda::getCudaEnabledDeviceCount() > 0)
829 seam_finder = makePtr<detail::GraphCutSeamFinderGpu>(GraphCutSeamFinderBase::COST_COLOR);
830 else
831#endif
832 seam_finder = makePtr<detail::GraphCutSeamFinder>(GraphCutSeamFinderBase::COST_COLOR);
833 }
834 else if (seam_find_type == "gc_colorgrad")
835 {
836#ifdef HAVE_OPENCV_CUDALEGACY
837 if (try_cuda && cuda::getCudaEnabledDeviceCount() > 0)
838 seam_finder = makePtr<detail::GraphCutSeamFinderGpu>(GraphCutSeamFinderBase::COST_COLOR_GRAD);
839 else
840#endif
841 seam_finder = makePtr<detail::GraphCutSeamFinder>(GraphCutSeamFinderBase::COST_COLOR_GRAD);
842 }
843 else if (seam_find_type == "dp_color")
844 seam_finder = makePtr<detail::DpSeamFinder>(DpSeamFinder::COLOR);
845 else if (seam_find_type == "dp_colorgrad")
846 seam_finder = makePtr<detail::DpSeamFinder>(DpSeamFinder::COLOR_GRAD);
847 if (!seam_finder)
848 {
849 cout << "Can't create the following seam finder '" << seam_find_type << "'\n";
850 return 1;
851 }
852
853 seam_finder->find(images_warped_f, corners, masks_warped);
854
855 LOGLN("Finding seams, time: " << ((getTickCount() - t) / getTickFrequency()) << " sec");
856
857 // Release unused memory
858 images.clear();
859 images_warped.clear();
860 images_warped_f.clear();
861 masks.clear();
862
863 LOGLN("Compositing...");
864#if ENABLE_LOG
865 t = getTickCount();
866#endif
867
868 Mat img_warped, img_warped_s;
869 Mat dilated_mask, seam_mask, mask, mask_warped;
870 Ptr<Blender> blender;
871 Ptr<Timelapser> timelapser;
872 //double compose_seam_aspect = 1;
873 double compose_work_aspect = 1;
874
875 for (int img_idx = 0; img_idx < num_images; ++img_idx)
876 {
877 LOGLN("Compositing image #" << indices[img_idx]+1);
878
879 // Read image and resize it if necessary
880 full_img = imread(samples::findFile(img_names[img_idx]));
881 if (!is_compose_scale_set)
882 {
883 if (compose_megapix > 0)
884 compose_scale = min(1.0, sqrt(compose_megapix * 1e6 / full_img.size().area()));
885 is_compose_scale_set = true;
886
887 // Compute relative scales
888 //compose_seam_aspect = compose_scale / seam_scale;
889 compose_work_aspect = compose_scale / work_scale;
890
891 // Update warped image scale
892 warped_image_scale *= static_cast<float>(compose_work_aspect);
893 warper = warper_creator->create(warped_image_scale);
894
895 // Update corners and sizes
896 for (int i = 0; i < num_images; ++i)
897 {
898 // Update intrinsics
899 cameras[i].focal *= compose_work_aspect;
900 cameras[i].ppx *= compose_work_aspect;
901 cameras[i].ppy *= compose_work_aspect;
902
903 // Update corner and size
904 Size sz = full_img_sizes[i];
905 if (std::abs(compose_scale - 1) > 1e-1)
906 {
907 sz.width = cvRound(full_img_sizes[i].width * compose_scale);
908 sz.height = cvRound(full_img_sizes[i].height * compose_scale);
909 }
910
911 Mat K;
912 cameras[i].K().convertTo(K, CV_32F);
913 Rect roi = warper->warpRoi(sz, K, cameras[i].R);
914 corners[i] = roi.tl();
915 sizes[i] = roi.size();
916 }
917 }
918 if (abs(compose_scale - 1) > 1e-1)
919 resize(full_img, img, Size(), compose_scale, compose_scale, INTER_LINEAR_EXACT);
920 else
921 img = full_img;
922 full_img.release();
923 Size img_size = img.size();
924
925 Mat K;
926 cameras[img_idx].K().convertTo(K, CV_32F);
927
928 // Warp the current image
929 warper->warp(img, K, cameras[img_idx].R, INTER_LINEAR, BORDER_REFLECT, img_warped);
930
931 // Warp the current image mask
932 mask.create(img_size, CV_8U);
933 mask.setTo(Scalar::all(255));
934 warper->warp(mask, K, cameras[img_idx].R, INTER_NEAREST, BORDER_CONSTANT, mask_warped);
935
936 // Compensate exposure
937 compensator->apply(img_idx, corners[img_idx], img_warped, mask_warped);
938
939 img_warped.convertTo(img_warped_s, CV_16S);
940 img_warped.release();
941 img.release();
942 mask.release();
943
944 dilate(masks_warped[img_idx], dilated_mask, Mat());
945 resize(dilated_mask, seam_mask, mask_warped.size(), 0, 0, INTER_LINEAR_EXACT);
946 mask_warped = seam_mask & mask_warped;
947
948 if (!blender && !timelapse)
949 {
950 blender = Blender::createDefault(blend_type, try_cuda);
951 Size dst_sz = resultRoi(corners, sizes).size();
952 float blend_width = sqrt(static_cast<float>(dst_sz.area())) * blend_strength / 100.f;
953 if (blend_width < 1.f)
954 blender = Blender::createDefault(Blender::NO, try_cuda);
955 else if (blend_type == Blender::MULTI_BAND)
956 {
957 MultiBandBlender* mb = dynamic_cast<MultiBandBlender*>(blender.get());
958 mb->setNumBands(static_cast<int>(ceil(log(blend_width)/log(2.)) - 1.));
959 LOGLN("Multi-band blender, number of bands: " << mb->numBands());
960 }
961 else if (blend_type == Blender::FEATHER)
962 {
963 FeatherBlender* fb = dynamic_cast<FeatherBlender*>(blender.get());
964 fb->setSharpness(1.f/blend_width);
965 LOGLN("Feather blender, sharpness: " << fb->sharpness());
966 }
967 blender->prepare(corners, sizes);
968 }
969 else if (!timelapser && timelapse)
970 {
971 timelapser = Timelapser::createDefault(timelapse_type);
972 timelapser->initialize(corners, sizes);
973 }
974
975 // Blend the current image
976 if (timelapse)
977 {
978 timelapser->process(img_warped_s, Mat::ones(img_warped_s.size(), CV_8UC1), corners[img_idx]);
979 String fixedFileName;
980 size_t pos_s = String(img_names[img_idx]).find_last_of("/\\");
981 if (pos_s == String::npos)
982 {
983 fixedFileName = "fixed_" + img_names[img_idx];
984 }
985 else
986 {
987 fixedFileName = "fixed_" + String(img_names[img_idx]).substr(pos_s + 1, String(img_names[img_idx]).length() - pos_s);
988 }
989 imwrite(fixedFileName, timelapser->getDst());
990 }
991 else
992 {
993 blender->feed(img_warped_s, mask_warped, corners[img_idx]);
994 }
995 }
996
997 if (!timelapse)
998 {
999 Mat result, result_mask;
1000 blender->blend(result, result_mask);
1001
1002 LOGLN("Compositing, time: " << ((getTickCount() - t) / getTickFrequency()) << " sec");
1003
1004 imwrite(result_name, result);
1005 }
1006
1007 LOGLN("Finished, total time: " << ((getTickCount() - app_start_time) / getTickFrequency()) << " sec");
1008 return 0;
1009}