1 /* Panorama_Tools - Generate, Edit and Convert Panoramic Images
2 Copyright (C) 1998,1999 - Helmut Dersch der@fh-furtwangen.de
4 This program is free software; you can redistribute it and/or modify
5 it under the terms of the GNU General Public License as published by
6 the Free Software Foundation; either version 2, or (at your option)
9 This program is distributed in the hope that it will be useful,
10 but WITHOUT ANY WARRANTY; without even the implied warranty of
11 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 GNU General Public License for more details.
14 You should have received a copy of the GNU General Public License
15 along with this program; if not, write to the Free Software
16 Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. */
18 /*------------------------------------------------------------*/
25 #include "utils_math.h"
26 #include "utils_image.h"
28 #define DEG_TO_RAD(x) ((x) * 2.0 * M_PI / 360.0 )
30 // Lookup Tables for Trig-functions and interpolator
35 static void matrix_matrix_mult( double m1[3][3],double m2[3][3],double result[3][3])
41 result[i][k] = m1[i][0] * m2[0][k] + m1[i][1] * m2[1][k] + m1[i][2] * m2[2][k];
42 } /* void matrix_matrix_mult */
44 /* Set matrix elements based on Euler angles a, b, c */
45 static void set_transformation_matrix (double m[3][3],
46 double a, double b, double c)
48 double mx[3][3], my[3][3], mz[3][3], dummy[3][3];
51 // Calculate Matrices;
53 mx[0][0] = 1.0 ; mx[0][1] = 0.0 ; mx[0][2] = 0.0;
54 mx[1][0] = 0.0 ; mx[1][1] = cos(a) ; mx[1][2] = sin(a);
55 mx[2][0] = 0.0 ; mx[2][1] =-mx[1][2] ; mx[2][2] = mx[1][1];
57 my[0][0] = cos(b); my[0][1] = 0.0 ; my[0][2] =-sin(b);
58 my[1][0] = 0.0 ; my[1][1] = 1.0 ; my[1][2] = 0.0;
59 my[2][0] = -my[0][2]; my[2][1] = 0.0 ; my[2][2] = my[0][0];
61 mz[0][0] = cos(c) ; mz[0][1] = sin(c) ; mz[0][2] = 0.0;
62 mz[1][0] =-mz[0][1] ; mz[1][1] = mz[0][0] ; mz[1][2] = 0.0;
63 mz[2][0] = 0.0 ; mz[2][1] = 0.0 ; mz[2][2] = 1.0;
65 /* Calculate `m = mz * mx * my' */
67 matrix_matrix_mult( mz, mx, dummy);
68 matrix_matrix_mult( dummy, my, m);
69 } /* void SetMatrix */
72 static double lanczos (double x)
74 static double table[257];
76 assert ((x >= 0.0) && (x <= 1.0));
83 for (i = 1; i <= 256; i++)
85 double t = ((double) i) / 256.0;
86 table[i] = 2 * sin (M_PI * t) * sin (M_PI_2 * t) / (M_PI * M_PI * t * t);
90 return (table[(int) (0.5 + 256.0 * x)]);
91 } /* double lanczos */
94 static int copy_pixel (ui_image_t *dest, const ui_image_t *src,
95 int x_dest, int y_dest,
96 double x_src_fp, double y_src_fp,
97 interpolator_t interp)
99 uint32_t pixel_dest = (y_dest * dest->width) + x_dest;
101 assert (x_src_fp >= 0.0);
106 if (interp == BILINEAR)
114 double y_bottom_frac;
118 x_src_left = (int) x_src_fp;
119 x_src_right = (x_src_left + 1) % src->width;
121 y_src_top = (int) y_src_fp;
122 y_src_bottom = y_src_top + 1;
123 if (y_src_bottom >= src->height)
124 y_src_bottom = src->height - 1;
126 x_right_frac = x_src_fp - x_src_left;
127 y_bottom_frac = y_src_fp - y_src_top;
130 * The polynomial function
132 * is similar to the Lanczos function, but much easier to compute. It's
133 * applied here to improve sharpness.
136 double t = x_src_fp - x_src_left;
137 x_right_frac = -2.0 * t * t * t + 3.0 * t * t;
139 t = y_src_fp - y_src_top;
140 y_bottom_frac = -2.0 * t * t * t + 3.0 * t * t;
143 assert ((x_right_frac >= 0.0) && (x_right_frac <= 1.0));
144 assert ((y_bottom_frac >= 0.0) && (y_bottom_frac <= 1.0));
146 for (i = 0; i < 3; i++)
148 uint8_t values[2][2];
153 values[0][0] = src->data[i][(y_src_top * src->width) + x_src_left];
154 values[0][1] = src->data[i][(y_src_top * src->width) + x_src_right];
155 values[1][0] = src->data[i][(y_src_bottom * src->width) + x_src_left];
156 values[1][1] = src->data[i][(y_src_bottom * src->width) + x_src_right];
158 value_left = (1.0 - y_bottom_frac) * values[0][0]
159 + y_bottom_frac * values[1][0];
160 value_right = (1.0 - y_bottom_frac) * values[0][1]
161 + y_bottom_frac * values[1][1];
163 value_final = (1.0 - x_right_frac) * value_left + x_right_frac * value_right;
165 assert ((value_final >= 0.0) && (value_final <= 255.0));
167 dest->data[i][pixel_dest] = (uint8_t) (value_final + 0.5);
170 else /* if (interp == NNEIGHBOUR) */
172 int x_src = (int) (x_src_fp + 0.5) % src->width;
173 int y_src = (int) (y_src_fp + 0.5) % src->height;
175 if ((x_src < 0) || (x_src >= src->width)
176 || (y_src < 0) || (y_src >= src->height))
178 dest->data[0][pixel_dest] = 0;
179 dest->data[1][pixel_dest] = 0;
180 dest->data[2][pixel_dest] = 0;
184 uint32_t pixel_src = (y_src * src->width) + x_src;
186 dest->data[0][pixel_dest] = src->data[0][pixel_src];
187 dest->data[1][pixel_dest] = src->data[1][pixel_src];
188 dest->data[2][pixel_dest] = src->data[2][pixel_src];
193 } /* int copy_pixel */
195 // Main transformation function. Destination image is calculated using transformation
196 // Function "func". Either all colors (color = 0) or one of rgb (color =1,2,3) are
197 // determined. If successful, TrPtr->success = 1. Memory for destination image
198 // must have been allocated and locked!
200 int pl_extract_view (ui_image_t *view, const ui_image_t *pano,
201 double pitch, double yaw, double fov, interpolator_t interp)
203 int x_dest, y_dest; // Loop through destination image
205 int dest_width_left = view->width / 2 ;
206 int dest_height_top = view->height / 2 ;
207 int src_width_left = pano->width / 2 ;
208 int src_height_top = pano->height / 2 ;
211 int x_min, x_max, y_min, y_max;
213 /* The transformation matrix */
219 { /* What the fuck does this? -octo */
220 double a = DEG_TO_RAD (fov);
221 double b = 2.0 * M_PI; /* DEG_TO_RAD (360.0) */
223 dist_r = ((double) view->width) / (2.0 * tan (a / 2.0));
224 dist_e = ((double) pano->width) / b;
227 set_transformation_matrix (tm, DEG_TO_RAD (pitch), DEG_TO_RAD (yaw), 0.0);
229 x_min = -dest_width_left; x_max = view->width - dest_width_left;
230 y_min = -dest_height_top; y_max = view->height - dest_height_top;
232 for(y_dest = y_min; y_dest < y_max; y_dest++)
234 for(x_dest = x_min; x_dest < x_max; x_dest++)
239 v[0] = tm[0][0] * x_dest + tm[1][0] * y_dest + tm[2][0] * dist_r;
240 v[1] = tm[0][1] * x_dest + tm[1][1] * y_dest + tm[2][1] * dist_r;
241 v[2] = tm[0][2] * x_dest + tm[1][2] * y_dest + tm[2][2] * dist_r;
243 x_src_fp = dist_e * atan2 (v[0], v[2]);
244 y_src_fp = dist_e * atan2 (v[1], sqrt (v[2] * v[2] + v[0] * v[0]));
246 copy_pixel (view, pano,
247 dest_width_left + x_dest, dest_height_top + y_dest,
248 src_width_left + x_src_fp, src_height_top + y_src_fp,
257 * vim: set tabstop=4 softtabstop=4 shiftwidth=4 :