1 //
2 // SPDX-License-Identifier: BSD-3-Clause
3 // Copyright (c) Contributors to the OpenEXR Project.
4 //
5 
6 
7 //-----------------------------------------------------------------------------
8 //
9 //	function blurImage() -- performs a hemispherical blur
10 //
11 //-----------------------------------------------------------------------------
12 
13 #include "blurImage.h"
14 
15 #include "namespaceAlias.h"
16 
17 #include <resizeImage.h>
18 #include <cstring>
19 #include "Iex.h"
20 #include <iostream>
21 #include <algorithm>
22 #include <string.h>
23 
24 
25 using namespace IMF;
26 using namespace std;
27 using namespace IMATH;
28 
29 
30 inline int
toInt(float x)31 toInt (float x)
32 {
33     return int (x + 0.5f);
34 }
35 
36 
37 inline double
sqr(double x)38 sqr (double x)
39 {
40     return x * x;
41 }
42 
43 
44 void
blurImage(EnvmapImage & image1,bool verbose)45 blurImage (EnvmapImage &image1, bool verbose)
46 {
47     //
48     // Ideally we would blur the input image directly by convolving
49     // it with a 180-degree wide blur kernel.  Unfortunately this
50     // is prohibitively expensive when the input image is large.
51     // In order to keep running times reasonable, we perform the
52     // blur on a small proxy image that will later be re-sampled
53     // to the desired output resolution.
54     //
55     // Here's how it works:
56     //
57     // * If the input image is in latitude-longitude format,
58     //   convert it into a cube-face environment map.
59     //
60     // * Repeatedly resample the image, each time shrinking
61     //   it to no less than half its current size, until the
62     //   width of each cube face is MAX_IN_WIDTH pixels.
63     //
64     // * Multiply each pixel by a weight that is proportinal
65     //   to the solid angle subtended by the pixel as seen
66     //   from the center of the environment cube.
67     //
68     // * Create an output image in cube-face format.
69     //   The cube faces of the output image are OUT_WIDTH
70     //   pixels wide.
71     //
72     // * For each pixel of the output image:
73     //
74     //       Set the output pixel's color to black
75     //
76     //       Determine the direction, d2, from the center of the
77     //       output environment cube to the center of the output
78     //	     pixel.
79     //
80     //       For each pixel of the input image:
81     //
82     //           Determine the direction, d1, from the center of
83     //           the input environment cube to the center of the
84     //           input pixel.
85     //
86     //           Multiply the input pixel's color by max (0, d1.dot(d2))
87     //           and add the result to the output pixel.
88     //
89 
90     const int MAX_IN_WIDTH = 40;
91     const int OUT_WIDTH = 100;
92 
93     if (verbose)
94 	cout << "blurring map image" << endl;
95 
96     EnvmapImage image2;
97     EnvmapImage *iptr1 = &image1;
98     EnvmapImage *iptr2 = &image2;
99 
100     int w = image1.dataWindow().max.x - image1.dataWindow().min.x + 1;
101     int h = w * 6;
102 
103     if (iptr1->type() == ENVMAP_LATLONG)
104     {
105 	//
106 	// Convert the input image from latitude-longitude
107 	// to cube-face format.
108 	//
109 
110 	if (verbose)
111 	    cout << "    converting to cube-face format" << endl;
112 
113 	w /= 4;
114 	h = w * 6;
115 
116 	Box2i dw (V2i (0, 0), V2i (w - 1, h - 1));
117 	resizeCube (*iptr1, *iptr2, dw, 1, 7);
118 
119 	swap (iptr1, iptr2);
120     }
121 
122     while (w > MAX_IN_WIDTH)
123     {
124 	//
125 	// Shrink the image.
126 	//
127 
128 	if (w >= MAX_IN_WIDTH * 2)
129 	    w /= 2;
130 	else
131 	    w = MAX_IN_WIDTH;
132 
133 	h = w * 6;
134 
135 	if (verbose)
136 	{
137 	    cout << "    resizing cube faces "
138 		    "to " << w << " by " << w << " pixels" << endl;
139 	}
140 
141 	Box2i dw (V2i (0, 0), V2i (w - 1, h - 1));
142 	resizeCube (*iptr1, *iptr2, dw, 1, 7);
143 
144 	swap (iptr1, iptr2);
145     }
146 
147     if (verbose)
148 	cout << "    computing pixel weights" << endl;
149 
150     {
151         //
152         // Multiply each pixel by a weight that is proportinal
153         // to the solid angle subtended by the pixel.
154         //
155 
156 	Box2i dw = iptr1->dataWindow();
157 	int sof = CubeMap::sizeOfFace (dw);
158 
159 	Array2D<Rgba> &pixels = iptr1->pixels();
160 
161 	double weightTotal = 0;
162 
163 	for (int f = CUBEFACE_POS_X; f <= CUBEFACE_NEG_Z; ++f)
164 	{
165 	    if (verbose)
166 		cout << "        face " << f << endl;
167 
168 	    CubeMapFace face = CubeMapFace (f);
169 	    V3f faceDir (0, 0, 0);
170             int ix = 0, iy = 0, iz = 0;
171 
172 	    switch (face)
173 	    {
174 	      case CUBEFACE_POS_X:
175 		faceDir = V3f (1, 0, 0);
176                 ix = 0;
177                 iy = 1;
178                 iz = 2;
179 		break;
180 
181 	      case CUBEFACE_NEG_X:
182 		faceDir = V3f (-1, 0, 0);
183                 ix = 0;
184                 iy = 1;
185                 iz = 2;
186 		break;
187 
188 	      case CUBEFACE_POS_Y:
189 		faceDir = V3f (0, 1, 0);
190                 ix = 1;
191                 iy = 0;
192                 iz = 2;
193 		break;
194 
195 	      case CUBEFACE_NEG_Y:
196 		faceDir = V3f (0, -1, 0);
197                 ix = 1;
198                 iy = 0;
199                 iz = 2;
200 		break;
201 
202 	      case CUBEFACE_POS_Z:
203 		faceDir = V3f (0, 0, 1);
204                 ix = 2;
205                 iy = 0;
206                 iz = 1;
207 		break;
208 
209 	      case CUBEFACE_NEG_Z:
210 		faceDir = V3f (0, 0, -1);
211                 ix = 2;
212                 iy = 0;
213                 iz = 1;
214 		break;
215 	    }
216 
217 	    for (int y = 0; y < sof; ++y)
218 	    {
219 		bool yEdge = (y == 0 || y == sof - 1);
220 
221 		for (int x = 0; x < sof; ++x)
222 		{
223 		    bool xEdge = (x == 0 || x == sof - 1);
224 
225 		    V2f posInFace (x, y);
226 
227 		    V3f dir =
228 			CubeMap::direction (face, dw, posInFace).normalized();
229 
230 		    V2f pos =
231 			CubeMap::pixelPosition (face, dw, posInFace);
232 
233                     //
234                     // The solid angle subtended by pixel (x,y), as seen
235                     // from the center of the cube, is proportional to the
236                     // square of the distance of the pixel from the center
237                     // of the cube and proportional to the dot product of
238                     // the viewing direction and the normal of the cube
239                     // face that contains the pixel.
240                     //
241 
242                     double weight =
243                         (dir ^ faceDir) *
244                         (sqr (dir[iy] / dir[ix]) + sqr (dir[iz] / dir[ix]) + 1);
245 
246                     //
247                     // Pixels at the edges and corners of the
248                     // cube are duplicated; we must adjust the
249                     // pixel weights accordingly.
250                     //
251 
252 		    if (xEdge && yEdge)
253 			weight /= 3;
254 		    else if (xEdge || yEdge)
255 			weight /= 2;
256 
257 		    Rgba &pixel = pixels[toInt (pos.y)][toInt (pos.x)];
258 
259 		    pixel.r *= weight;
260 		    pixel.g *= weight;
261 		    pixel.b *= weight;
262 		    pixel.a *= weight;
263 
264 		    weightTotal += weight;
265 		}
266 	    }
267 	}
268 
269 	//
270 	// The weighting operation above has made the overall image darker.
271 	// Apply a correction to recover the image's original brightness.
272 	//
273 
274 	int w = dw.max.x - dw.min.x + 1;
275 	int h = dw.max.y - dw.min.y + 1;
276 	size_t numPixels = w * h;
277 	double weight = numPixels / weightTotal;
278 
279 	Rgba *p = &pixels[0][0];
280 	Rgba *end = p + numPixels;
281 
282 	while (p < end)
283 	{
284 	    p->r *= weight;
285 	    p->g *= weight;
286 	    p->b *= weight;
287 	    p->a *= weight;
288 
289 	    ++p;
290 	}
291     }
292 
293     {
294 	if (verbose)
295 	    cout << "    generating blurred image" << endl;
296 
297 	Box2i dw1 = iptr1->dataWindow();
298 	int sof1 = CubeMap::sizeOfFace (dw1);
299 
300 	Box2i dw2 (V2i (0, 0), V2i (OUT_WIDTH - 1, OUT_WIDTH * 6 - 1));
301 	int sof2 = CubeMap::sizeOfFace (dw2);
302 
303 	iptr2->resize (ENVMAP_CUBE, dw2);
304 	iptr2->clear ();
305 
306 	Array2D<Rgba> &pixels1 = iptr1->pixels();
307 	Array2D<Rgba> &pixels2 = iptr2->pixels();
308 
309 	for (int f2 = CUBEFACE_POS_X; f2 <= CUBEFACE_NEG_Z; ++f2)
310 	{
311 	    if (verbose)
312 		cout << "        face " << f2 << endl;
313 
314 	    CubeMapFace face2 = CubeMapFace (f2);
315 
316 	    for (int y2 = 0; y2 < sof2; ++y2)
317 	    {
318 		for (int x2 = 0; x2 < sof2; ++x2)
319 		{
320 		    V2f posInFace2 (x2, y2);
321 
322 		    V3f dir2 = CubeMap::direction
323 			(face2, dw2, posInFace2);
324 
325 		    V2f pos2 = CubeMap::pixelPosition
326 			(face2, dw2, posInFace2);
327 
328 		    double weightTotal = 0;
329 		    double rTotal = 0;
330 		    double gTotal = 0;
331 		    double bTotal = 0;
332 		    double aTotal = 0;
333 
334 		    Rgba &pixel2 =
335 			pixels2[toInt (pos2.y)][toInt (pos2.x)];
336 
337 		    for (int f1 = CUBEFACE_POS_X; f1 <= CUBEFACE_NEG_Z; ++f1)
338 		    {
339 			CubeMapFace face1 = CubeMapFace (f1);
340 
341 			for (int y1 = 0; y1 < sof1; ++y1)
342 			{
343 			    for (int x1 = 0; x1 < sof1; ++x1)
344 			    {
345 				V2f posInFace1 (x1, y1);
346 
347 				V3f dir1 = CubeMap::direction
348 				    (face1, dw1, posInFace1);
349 
350 				V2f pos1 = CubeMap::pixelPosition
351 				    (face1, dw1, posInFace1);
352 
353 				double weight = dir1 ^ dir2;
354 
355 				if (weight <= 0)
356 				    continue;
357 
358 				Rgba &pixel1 =
359 				    pixels1[toInt (pos1.y)][toInt (pos1.x)];
360 
361 				weightTotal += weight;
362 				rTotal += pixel1.r * weight;
363 				gTotal += pixel1.g * weight;
364 				bTotal += pixel1.b * weight;
365 				aTotal += pixel1.a * weight;
366 			    }
367 			}
368 		    }
369 
370 		    pixel2.r = rTotal / weightTotal;
371 		    pixel2.g = gTotal / weightTotal;
372 		    pixel2.b = bTotal / weightTotal;
373 		    pixel2.a = aTotal / weightTotal;
374 		}
375 	    }
376 	}
377 
378 	swap (iptr1, iptr2);
379     }
380 
381     //
382     // Depending on how many times we've re-sampled the image,
383     // the result is now either in image1 or in image2.
384     // If necessary, copy the result into image1.
385     //
386 
387     if (iptr1 != &image1)
388     {
389 	if (verbose)
390 	    cout << "    copying" << endl;
391 
392 	Box2i dw = iptr1->dataWindow();
393 	image1.resize (ENVMAP_CUBE, dw);
394 
395 	int w = dw.max.x - dw.min.x + 1;
396 	int h = dw.max.y - dw.min.y + 1;
397 	size_t size = w * h * sizeof (Rgba);
398 
399 	memcpy (&image1.pixels()[0][0], &iptr1->pixels()[0][0], size);
400     }
401 }
402