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 #include "opencl_kernels_stitching.hpp"
45 
46 namespace cv {
47 namespace detail {
48 
setCameraParams(InputArray _K,InputArray _R,InputArray _T)49 void ProjectorBase::setCameraParams(InputArray _K, InputArray _R, InputArray _T)
50 {
51     Mat K = _K.getMat(), R = _R.getMat(), T = _T.getMat();
52 
53     CV_Assert(K.size() == Size(3, 3) && K.type() == CV_32F);
54     CV_Assert(R.size() == Size(3, 3) && R.type() == CV_32F);
55     CV_Assert((T.size() == Size(1, 3) || T.size() == Size(3, 1)) && T.type() == CV_32F);
56 
57     Mat_<float> K_(K);
58     k[0] = K_(0,0); k[1] = K_(0,1); k[2] = K_(0,2);
59     k[3] = K_(1,0); k[4] = K_(1,1); k[5] = K_(1,2);
60     k[6] = K_(2,0); k[7] = K_(2,1); k[8] = K_(2,2);
61 
62     Mat_<float> Rinv = R.t();
63     rinv[0] = Rinv(0,0); rinv[1] = Rinv(0,1); rinv[2] = Rinv(0,2);
64     rinv[3] = Rinv(1,0); rinv[4] = Rinv(1,1); rinv[5] = Rinv(1,2);
65     rinv[6] = Rinv(2,0); rinv[7] = Rinv(2,1); rinv[8] = Rinv(2,2);
66 
67     Mat_<float> R_Kinv = R * K.inv();
68     r_kinv[0] = R_Kinv(0,0); r_kinv[1] = R_Kinv(0,1); r_kinv[2] = R_Kinv(0,2);
69     r_kinv[3] = R_Kinv(1,0); r_kinv[4] = R_Kinv(1,1); r_kinv[5] = R_Kinv(1,2);
70     r_kinv[6] = R_Kinv(2,0); r_kinv[7] = R_Kinv(2,1); r_kinv[8] = R_Kinv(2,2);
71 
72     Mat_<float> K_Rinv = K * Rinv;
73     k_rinv[0] = K_Rinv(0,0); k_rinv[1] = K_Rinv(0,1); k_rinv[2] = K_Rinv(0,2);
74     k_rinv[3] = K_Rinv(1,0); k_rinv[4] = K_Rinv(1,1); k_rinv[5] = K_Rinv(1,2);
75     k_rinv[6] = K_Rinv(2,0); k_rinv[7] = K_Rinv(2,1); k_rinv[8] = K_Rinv(2,2);
76 
77     Mat_<float> T_(T.reshape(0, 3));
78     t[0] = T_(0,0); t[1] = T_(1,0); t[2] = T_(2,0);
79 }
80 
81 
warpPoint(const Point2f & pt,InputArray K,InputArray R,InputArray T)82 Point2f PlaneWarper::warpPoint(const Point2f &pt, InputArray K, InputArray R, InputArray T)
83 {
84     projector_.setCameraParams(K, R, T);
85     Point2f uv;
86     projector_.mapForward(pt.x, pt.y, uv.x, uv.y);
87     return uv;
88 }
89 
warpPoint(const Point2f & pt,InputArray K,InputArray R)90 Point2f PlaneWarper::warpPoint(const Point2f &pt, InputArray K, InputArray R)
91 {
92     float tz[] = {0.f, 0.f, 0.f};
93     Mat_<float> T(3, 1, tz);
94     return warpPoint(pt, K, R, T);
95 }
96 
buildMaps(Size src_size,InputArray K,InputArray R,OutputArray xmap,OutputArray ymap)97 Rect PlaneWarper::buildMaps(Size src_size, InputArray K, InputArray R, OutputArray xmap, OutputArray ymap)
98 {
99     return buildMaps(src_size, K, R, Mat::zeros(3, 1, CV_32FC1), xmap, ymap);
100 }
101 
buildMaps(Size src_size,InputArray K,InputArray R,InputArray T,OutputArray _xmap,OutputArray _ymap)102 Rect PlaneWarper::buildMaps(Size src_size, InputArray K, InputArray R, InputArray T, OutputArray _xmap, OutputArray _ymap)
103 {
104     projector_.setCameraParams(K, R, T);
105 
106     Point dst_tl, dst_br;
107     detectResultRoi(src_size, dst_tl, dst_br);
108 
109     Size dsize(dst_br.x - dst_tl.x + 1, dst_br.y - dst_tl.y + 1);
110     _xmap.create(dsize, CV_32FC1);
111     _ymap.create(dsize, CV_32FC1);
112 
113     if (ocl::useOpenCL())
114     {
115         ocl::Kernel k("buildWarpPlaneMaps", ocl::stitching::warpers_oclsrc);
116         if (!k.empty())
117         {
118             int rowsPerWI = ocl::Device::getDefault().isIntel() ? 4 : 1;
119             Mat k_rinv(1, 9, CV_32FC1, projector_.k_rinv), t(1, 3, CV_32FC1, projector_.t);
120             UMat uxmap = _xmap.getUMat(), uymap = _ymap.getUMat(),
121                     uk_rinv = k_rinv.getUMat(ACCESS_READ), ut = t.getUMat(ACCESS_READ);
122 
123             k.args(ocl::KernelArg::WriteOnlyNoSize(uxmap), ocl::KernelArg::WriteOnly(uymap),
124                    ocl::KernelArg::PtrReadOnly(uk_rinv), ocl::KernelArg::PtrReadOnly(ut),
125                    dst_tl.x, dst_tl.y, 1/projector_.scale, rowsPerWI);
126 
127             size_t globalsize[2] = { dsize.width, (dsize.height + rowsPerWI - 1) / rowsPerWI };
128             if (k.run(2, globalsize, NULL, true))
129             {
130                 CV_IMPL_ADD(CV_IMPL_OCL);
131                 return Rect(dst_tl, dst_br);
132             }
133         }
134     }
135 
136     Mat xmap = _xmap.getMat(), ymap = _ymap.getMat();
137 
138     float x, y;
139     for (int v = dst_tl.y; v <= dst_br.y; ++v)
140     {
141         for (int u = dst_tl.x; u <= dst_br.x; ++u)
142         {
143             projector_.mapBackward(static_cast<float>(u), static_cast<float>(v), x, y);
144             xmap.at<float>(v - dst_tl.y, u - dst_tl.x) = x;
145             ymap.at<float>(v - dst_tl.y, u - dst_tl.x) = y;
146         }
147     }
148 
149     return Rect(dst_tl, dst_br);
150 }
151 
152 
warp(InputArray src,InputArray K,InputArray R,InputArray T,int interp_mode,int border_mode,OutputArray dst)153 Point PlaneWarper::warp(InputArray src, InputArray K, InputArray R, InputArray T, int interp_mode, int border_mode,
154                         OutputArray dst)
155 {
156     UMat uxmap, uymap;
157     Rect dst_roi = buildMaps(src.size(), K, R, T, uxmap, uymap);
158 
159     dst.create(dst_roi.height + 1, dst_roi.width + 1, src.type());
160     remap(src, dst, uxmap, uymap, interp_mode, border_mode);
161 
162     return dst_roi.tl();
163 }
164 
warp(InputArray src,InputArray K,InputArray R,int interp_mode,int border_mode,OutputArray dst)165 Point PlaneWarper::warp(InputArray src, InputArray K, InputArray R,
166                         int interp_mode, int border_mode, OutputArray dst)
167 {
168     float tz[] = {0.f, 0.f, 0.f};
169     Mat_<float> T(3, 1, tz);
170     return warp(src, K, R, T, interp_mode, border_mode, dst);
171 }
172 
warpRoi(Size src_size,InputArray K,InputArray R,InputArray T)173 Rect PlaneWarper::warpRoi(Size src_size, InputArray K, InputArray R, InputArray T)
174 {
175     projector_.setCameraParams(K, R, T);
176 
177     Point dst_tl, dst_br;
178     detectResultRoi(src_size, dst_tl, dst_br);
179 
180     return Rect(dst_tl, Point(dst_br.x + 1, dst_br.y + 1));
181 }
182 
warpRoi(Size src_size,InputArray K,InputArray R)183 Rect PlaneWarper::warpRoi(Size src_size, InputArray K, InputArray R)
184 {
185     float tz[] = {0.f, 0.f, 0.f};
186     Mat_<float> T(3, 1, tz);
187     return warpRoi(src_size, K, R, T);
188 }
189 
190 
detectResultRoi(Size src_size,Point & dst_tl,Point & dst_br)191 void PlaneWarper::detectResultRoi(Size src_size, Point &dst_tl, Point &dst_br)
192 {
193     float tl_uf = std::numeric_limits<float>::max();
194     float tl_vf = std::numeric_limits<float>::max();
195     float br_uf = -std::numeric_limits<float>::max();
196     float br_vf = -std::numeric_limits<float>::max();
197 
198     float u, v;
199 
200     projector_.mapForward(0, 0, u, v);
201     tl_uf = std::min(tl_uf, u); tl_vf = std::min(tl_vf, v);
202     br_uf = std::max(br_uf, u); br_vf = std::max(br_vf, v);
203 
204     projector_.mapForward(0, static_cast<float>(src_size.height - 1), u, v);
205     tl_uf = std::min(tl_uf, u); tl_vf = std::min(tl_vf, v);
206     br_uf = std::max(br_uf, u); br_vf = std::max(br_vf, v);
207 
208     projector_.mapForward(static_cast<float>(src_size.width - 1), 0, u, v);
209     tl_uf = std::min(tl_uf, u); tl_vf = std::min(tl_vf, v);
210     br_uf = std::max(br_uf, u); br_vf = std::max(br_vf, v);
211 
212     projector_.mapForward(static_cast<float>(src_size.width - 1), static_cast<float>(src_size.height - 1), u, v);
213     tl_uf = std::min(tl_uf, u); tl_vf = std::min(tl_vf, v);
214     br_uf = std::max(br_uf, u); br_vf = std::max(br_vf, v);
215 
216     dst_tl.x = static_cast<int>(tl_uf);
217     dst_tl.y = static_cast<int>(tl_vf);
218     dst_br.x = static_cast<int>(br_uf);
219     dst_br.y = static_cast<int>(br_vf);
220 }
221 
222 
detectResultRoi(Size src_size,Point & dst_tl,Point & dst_br)223 void SphericalWarper::detectResultRoi(Size src_size, Point &dst_tl, Point &dst_br)
224 {
225     detectResultRoiByBorder(src_size, dst_tl, dst_br);
226 
227     float tl_uf = static_cast<float>(dst_tl.x);
228     float tl_vf = static_cast<float>(dst_tl.y);
229     float br_uf = static_cast<float>(dst_br.x);
230     float br_vf = static_cast<float>(dst_br.y);
231 
232     float x = projector_.rinv[1];
233     float y = projector_.rinv[4];
234     float z = projector_.rinv[7];
235     if (y > 0.f)
236     {
237         float x_ = (projector_.k[0] * x + projector_.k[1] * y) / z + projector_.k[2];
238         float y_ = projector_.k[4] * y / z + projector_.k[5];
239         if (x_ > 0.f && x_ < src_size.width && y_ > 0.f && y_ < src_size.height)
240         {
241             tl_uf = std::min(tl_uf, 0.f); tl_vf = std::min(tl_vf, static_cast<float>(CV_PI * projector_.scale));
242             br_uf = std::max(br_uf, 0.f); br_vf = std::max(br_vf, static_cast<float>(CV_PI * projector_.scale));
243         }
244     }
245 
246     x = projector_.rinv[1];
247     y = -projector_.rinv[4];
248     z = projector_.rinv[7];
249     if (y > 0.f)
250     {
251         float x_ = (projector_.k[0] * x + projector_.k[1] * y) / z + projector_.k[2];
252         float y_ = projector_.k[4] * y / z + projector_.k[5];
253         if (x_ > 0.f && x_ < src_size.width && y_ > 0.f && y_ < src_size.height)
254         {
255             tl_uf = std::min(tl_uf, 0.f); tl_vf = std::min(tl_vf, static_cast<float>(0));
256             br_uf = std::max(br_uf, 0.f); br_vf = std::max(br_vf, static_cast<float>(0));
257         }
258     }
259 
260     dst_tl.x = static_cast<int>(tl_uf);
261     dst_tl.y = static_cast<int>(tl_vf);
262     dst_br.x = static_cast<int>(br_uf);
263     dst_br.y = static_cast<int>(br_vf);
264 }
265 
detectResultRoi(Size src_size,Point & dst_tl,Point & dst_br)266 void SphericalPortraitWarper::detectResultRoi(Size src_size, Point &dst_tl, Point &dst_br)
267 {
268     detectResultRoiByBorder(src_size, dst_tl, dst_br);
269 
270     float tl_uf = static_cast<float>(dst_tl.x);
271     float tl_vf = static_cast<float>(dst_tl.y);
272     float br_uf = static_cast<float>(dst_br.x);
273     float br_vf = static_cast<float>(dst_br.y);
274 
275     float x = projector_.rinv[0];
276     float y = projector_.rinv[3];
277     float z = projector_.rinv[6];
278     if (y > 0.f)
279     {
280         float x_ = (projector_.k[0] * x + projector_.k[1] * y) / z + projector_.k[2];
281         float y_ = projector_.k[4] * y / z + projector_.k[5];
282         if (x_ > 0.f && x_ < src_size.width && y_ > 0.f && y_ < src_size.height)
283         {
284             tl_uf = std::min(tl_uf, 0.f); tl_vf = std::min(tl_vf, static_cast<float>(CV_PI * projector_.scale));
285             br_uf = std::max(br_uf, 0.f); br_vf = std::max(br_vf, static_cast<float>(CV_PI * projector_.scale));
286         }
287     }
288 
289     x = projector_.rinv[0];
290     y = -projector_.rinv[3];
291     z = projector_.rinv[6];
292     if (y > 0.f)
293     {
294         float x_ = (projector_.k[0] * x + projector_.k[1] * y) / z + projector_.k[2];
295         float y_ = projector_.k[4] * y / z + projector_.k[5];
296         if (x_ > 0.f && x_ < src_size.width && y_ > 0.f && y_ < src_size.height)
297         {
298             tl_uf = std::min(tl_uf, 0.f); tl_vf = std::min(tl_vf, static_cast<float>(0));
299             br_uf = std::max(br_uf, 0.f); br_vf = std::max(br_vf, static_cast<float>(0));
300         }
301     }
302 
303     dst_tl.x = static_cast<int>(tl_uf);
304     dst_tl.y = static_cast<int>(tl_vf);
305     dst_br.x = static_cast<int>(br_uf);
306     dst_br.y = static_cast<int>(br_vf);
307 }
308 
309 /////////////////////////////////////////// SphericalWarper ////////////////////////////////////////
310 
buildMaps(Size src_size,InputArray K,InputArray R,OutputArray xmap,OutputArray ymap)311 Rect SphericalWarper::buildMaps(Size src_size, InputArray K, InputArray R, OutputArray xmap, OutputArray ymap)
312 {
313     if (ocl::useOpenCL())
314     {
315         ocl::Kernel k("buildWarpSphericalMaps", ocl::stitching::warpers_oclsrc);
316         if (!k.empty())
317         {
318             int rowsPerWI = ocl::Device::getDefault().isIntel() ? 4 : 1;
319             projector_.setCameraParams(K, R);
320 
321             Point dst_tl, dst_br;
322             detectResultRoi(src_size, dst_tl, dst_br);
323 
324             Size dsize(dst_br.x - dst_tl.x + 1, dst_br.y - dst_tl.y + 1);
325             xmap.create(dsize, CV_32FC1);
326             ymap.create(dsize, CV_32FC1);
327 
328             Mat k_rinv(1, 9, CV_32FC1, projector_.k_rinv);
329             UMat uxmap = xmap.getUMat(), uymap = ymap.getUMat(), uk_rinv = k_rinv.getUMat(ACCESS_READ);
330 
331             k.args(ocl::KernelArg::WriteOnlyNoSize(uxmap), ocl::KernelArg::WriteOnly(uymap),
332                    ocl::KernelArg::PtrReadOnly(uk_rinv), dst_tl.x, dst_tl.y, 1/projector_.scale, rowsPerWI);
333 
334             size_t globalsize[2] = { dsize.width, (dsize.height + rowsPerWI - 1) / rowsPerWI };
335             if (k.run(2, globalsize, NULL, true))
336             {
337                 CV_IMPL_ADD(CV_IMPL_OCL);
338                 return Rect(dst_tl, dst_br);
339             }
340         }
341     }
342 
343     return RotationWarperBase<SphericalProjector>::buildMaps(src_size, K, R, xmap, ymap);
344 }
345 
warp(InputArray src,InputArray K,InputArray R,int interp_mode,int border_mode,OutputArray dst)346 Point SphericalWarper::warp(InputArray src, InputArray K, InputArray R, int interp_mode, int border_mode, OutputArray dst)
347 {
348     UMat uxmap, uymap;
349     Rect dst_roi = buildMaps(src.size(), K, R, uxmap, uymap);
350 
351     dst.create(dst_roi.height + 1, dst_roi.width + 1, src.type());
352     remap(src, dst, uxmap, uymap, interp_mode, border_mode);
353 
354     return dst_roi.tl();
355 }
356 
357 /////////////////////////////////////////// CylindricalWarper ////////////////////////////////////////
358 
buildMaps(Size src_size,InputArray K,InputArray R,OutputArray xmap,OutputArray ymap)359 Rect CylindricalWarper::buildMaps(Size src_size, InputArray K, InputArray R, OutputArray xmap, OutputArray ymap)
360 {
361     if (ocl::useOpenCL())
362     {
363         ocl::Kernel k("buildWarpCylindricalMaps", ocl::stitching::warpers_oclsrc);
364         if (!k.empty())
365         {
366             int rowsPerWI = ocl::Device::getDefault().isIntel() ? 4 : 1;
367             projector_.setCameraParams(K, R);
368 
369             Point dst_tl, dst_br;
370             detectResultRoi(src_size, dst_tl, dst_br);
371 
372             Size dsize(dst_br.x - dst_tl.x + 1, dst_br.y - dst_tl.y + 1);
373             xmap.create(dsize, CV_32FC1);
374             ymap.create(dsize, CV_32FC1);
375 
376             Mat k_rinv(1, 9, CV_32FC1, projector_.k_rinv);
377             UMat uxmap = xmap.getUMat(), uymap = ymap.getUMat(), uk_rinv = k_rinv.getUMat(ACCESS_READ);
378 
379             k.args(ocl::KernelArg::WriteOnlyNoSize(uxmap), ocl::KernelArg::WriteOnly(uymap),
380                    ocl::KernelArg::PtrReadOnly(uk_rinv), dst_tl.x, dst_tl.y, 1/projector_.scale,
381                    rowsPerWI);
382 
383             size_t globalsize[2] = { dsize.width, (dsize.height + rowsPerWI - 1) / rowsPerWI };
384             if (k.run(2, globalsize, NULL, true))
385             {
386                 CV_IMPL_ADD(CV_IMPL_OCL);
387                 return Rect(dst_tl, dst_br);
388             }
389         }
390     }
391 
392     return RotationWarperBase<CylindricalProjector>::buildMaps(src_size, K, R, xmap, ymap);
393 }
394 
warp(InputArray src,InputArray K,InputArray R,int interp_mode,int border_mode,OutputArray dst)395 Point CylindricalWarper::warp(InputArray src, InputArray K, InputArray R, int interp_mode, int border_mode, OutputArray dst)
396 {
397     UMat uxmap, uymap;
398     Rect dst_roi = buildMaps(src.size(), K, R, uxmap, uymap);
399 
400     dst.create(dst_roi.height + 1, dst_roi.width + 1, src.type());
401     remap(src, dst, uxmap, uymap, interp_mode, border_mode);
402 
403     return dst_roi.tl();
404 }
405 
406 } // namespace detail
407 } // namespace cv
408