1 /* -*- c -*- */
2
3 /*
4 * zoom.c
5 *
6 * metapixel
7 *
8 * Copyright (C) 2004 Mark Probst
9 *
10 * This program is free software; you can redistribute it and/or
11 * modify it under the terms of the GNU General Public License
12 * as published by the Free Software Foundation; either version 2
13 * of the License, or (at your option) any later version.
14 *
15 * This program is distributed in the hope that it will be useful,
16 * but WITHOUT ANY WARRANTY; without even the implied warranty of
17 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 * GNU General Public License for more details.
19 *
20 * You should have received a copy of the GNU General Public License
21 * along with this program; if not, write to the Free Software
22 * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
23 */
24
25 #include <math.h>
26 #include <stdlib.h>
27 #include <assert.h>
28
29 #include "zoom.h"
30
31 #ifndef MIN
32 #define MIN(a,b) ((a)<(b)?(a):(b))
33 #endif
34 #ifndef MAX
35 #define MAX(a,b) ((a)>(b)?(a):(b))
36 #endif
37
38 #define MAX_FILTER FILTER_MITCHELL
39
40 typedef struct
41 {
42 int index;
43 float weight;
44 int iweight;
45 } sample_t;
46
47 typedef struct
48 {
49 int num_samples;
50 sample_t samples[0];
51 } sample_window_t;
52
53 static float
filter_box(float x)54 filter_box (float x)
55 {
56 if (x < -0.5)
57 return 0.0;
58 if (x <= 0.5)
59 return 1.0;
60 return 0.0;
61 }
62
63 static float
filter_triangle(float x)64 filter_triangle (float x)
65 {
66 if (x < -1.0)
67 return 0.0;
68 if (x < 0.0)
69 return 1.0 + x;
70 if (x < 1.0)
71 return 1.0 - x;
72 return 0.0;
73 }
74
75 /*
76 * see Mitchell&Netravali, "Reconstruction Filters in Computer
77 * Graphics", Proceedings of the 15th annual conference on Computer
78 * graphics and interactive techniques, ACM Press 1988
79 */
80 static float
filter_mitchell(float x)81 filter_mitchell (float x)
82 {
83 #define B (1.0 / 3.0)
84 #define C (1.0 / 3.0)
85
86 static float a0 = ( 6.0 - 2.0 * B ) / 6.0;
87 static float a2 = (-18.0 + 12.0 * B + 6.0 * C) / 6.0;
88 static float a3 = ( 12.0 - 9.0 * B - 6.0 * C) / 6.0;
89
90 static float b0 = ( 8.0 * B + 24.0 * C) / 6.0;
91 static float b1 = ( - 12.0 * B - 48.0 * C) / 6.0;
92 static float b2 = ( 6.0 * B + 30.0 * C) / 6.0;
93 static float b3 = ( - B - 6.0 * C) / 6.0;
94
95 x = fabsf(x);
96
97 if (x < 1.0)
98 return a0 + (x * x) * (a2 + x * a3);
99 if (x < 2.0)
100 return b0 + x * (b1 + x * (b2 + x * b3));
101 return 0.0;
102
103 #undef b
104 #undef c
105 }
106
107 static filter_t filters[] = {
108 { &filter_box, 0.5 },
109 { &filter_triangle, 1.0 },
110 { &filter_mitchell, 2.0 }
111 };
112
113 filter_t*
get_filter(int index)114 get_filter (int index)
115 {
116 if (index < 0 || index > MAX_FILTER)
117 return 0;
118
119 return &filters[index];
120 }
121
122 #define NUM_ACCURACY_BITS 12
123
124 static sample_window_t*
make_sample_window(float center,float scale,float support_radius,filter_func_t filter_func,int num_indexes)125 make_sample_window (float center, float scale, float support_radius, filter_func_t filter_func, int num_indexes)
126 {
127 float lower_bound = center - support_radius;
128 float upper_bound = center + support_radius;
129 int lower_index = floor(lower_bound + 0.5);
130 int upper_index = floor(upper_bound - 0.5);
131 int num_samples;
132 sample_window_t *window;
133 int i;
134 float weight_sum;
135
136 lower_index = MAX(0, lower_index);
137 upper_index = MIN(num_indexes - 1, upper_index);
138
139 if (upper_index < lower_index)
140 upper_index = lower_index = floor(center);
141
142 num_samples = upper_index - lower_index + 1;
143
144 assert(num_samples > 0);
145
146 window = (sample_window_t*)malloc(sizeof(sample_window_t) + num_samples * sizeof(sample_t));
147 assert(window != 0);
148
149 window->num_samples = num_samples;
150
151 weight_sum = 0.0;
152 for (i = 0; i < num_samples; ++i)
153 {
154 int index = lower_index + i;
155 float sample_center = (float)index + 0.5;
156
157 window->samples[i].index = index;
158 window->samples[i].weight = filter_func((sample_center - center) / scale);
159
160 weight_sum += window->samples[i].weight;
161 }
162
163 assert(weight_sum > 0.0);
164
165 for (i = 0; i < num_samples; ++i)
166 {
167 window->samples[i].weight /= weight_sum;
168 window->samples[i].iweight = (1 << NUM_ACCURACY_BITS) * window->samples[i].weight;
169 }
170
171 return window;
172 }
173
174 static sample_window_t**
make_sample_windows(float filter_scale,float filter_support_radius,filter_func_t filter_func,int dest_size,int src_size,float scale)175 make_sample_windows (float filter_scale, float filter_support_radius, filter_func_t filter_func,
176 int dest_size, int src_size, float scale)
177 {
178 sample_window_t **sample_windows;
179 int i;
180
181 sample_windows = (sample_window_t**)malloc(dest_size * sizeof(sample_window_t*));
182 assert(sample_windows != 0);
183 for (i = 0; i < dest_size; ++i)
184 {
185 float dest_center = (float)i + 0.5;
186 float src_center = dest_center / scale;
187
188 sample_windows[i] = make_sample_window(src_center, filter_scale, filter_support_radius,
189 filter_func, src_size);
190 assert(sample_windows[i] != 0);
191 }
192
193 return sample_windows;
194 }
195
196 static void
free_sample_windows(sample_window_t ** sample_windows,int size)197 free_sample_windows (sample_window_t **sample_windows, int size)
198 {
199 int i;
200
201 for (i = 0; i < size; ++i)
202 free(sample_windows[i]);
203 free(sample_windows);
204 }
205
206 static void
zoom_unidirectional(unsigned char * dest,unsigned char * src,int num_channels,sample_window_t ** sample_windows,int num_pixels_in_entity,int num_entities,int dest_pixel_advance,int src_pixel_advance,int dest_entity_advance,int src_entity_advance)207 zoom_unidirectional (unsigned char *dest, unsigned char *src, int num_channels, sample_window_t **sample_windows,
208 int num_pixels_in_entity, int num_entities,
209 int dest_pixel_advance, int src_pixel_advance,
210 int dest_entity_advance, int src_entity_advance)
211 {
212 int i;
213 unsigned char *dest_entity, *src_entity;
214 int channels[num_channels];
215
216 dest_entity = dest;
217 src_entity = src;
218 for (i = 0; i < num_entities; ++i)
219 {
220 int j;
221 unsigned char *dest_pixel;
222
223 dest_pixel = dest_entity;
224 for (j = 0; j < num_pixels_in_entity; ++j)
225 {
226
227 int k;
228
229 for (k = 0; k < num_channels; ++k)
230 channels[k] = 0;
231
232 for (k = 0; k < sample_windows[j]->num_samples; ++k)
233 {
234 int l;
235 sample_t *sample = &sample_windows[j]->samples[k];
236 unsigned char *src_pixel = &src_entity[sample->index * src_pixel_advance];
237
238 for (l = 0; l < num_channels; ++l)
239 channels[l] += (int)src_pixel[l] * sample->iweight;
240 /* ((((int)src_pixel[l]) << CHANNEL_SHIFT) + (1 << (CHANNEL_SHIFT - 1))) * sample->weight; */
241 }
242
243 for (k = 0; k < num_channels; ++k)
244 {
245 int value = channels[k] >> NUM_ACCURACY_BITS;
246
247 dest_pixel[k] = MAX(0, MIN(255, value));
248 }
249
250 dest_pixel += dest_pixel_advance;
251 }
252
253 dest_entity += dest_entity_advance;
254 src_entity += src_entity_advance;
255 }
256 }
257
258 void
zoom_image(unsigned char * dest,unsigned char * src,filter_t * filter,int num_channels,int dest_width,int dest_height,int dest_row_stride,int src_width,int src_height,int src_row_stride)259 zoom_image (unsigned char *dest, unsigned char *src,
260 filter_t *filter, int num_channels,
261 int dest_width, int dest_height, int dest_row_stride,
262 int src_width, int src_height, int src_row_stride)
263 {
264 float x_scale, y_scale;
265 float filter_x_scale, filter_y_scale;
266 float filter_x_support_radius, filter_y_support_radius;
267 sample_window_t **x_sample_windows, **y_sample_windows;
268 unsigned char *temp_image;
269
270 assert(dest != 0 && src != 0 && filter != 0);
271
272 assert(dest_width > 0 && dest_height > 0);
273
274 x_scale = (float)dest_width / (float)src_width;
275 y_scale = (float)dest_height / (float)src_height;
276
277 filter_x_scale = MAX(1.0, 1.0 / x_scale);
278 filter_y_scale = MAX(1.0, 1.0 / y_scale);
279
280 filter_x_support_radius = filter->support_radius * filter_x_scale;
281 filter_y_support_radius = filter->support_radius * filter_y_scale;
282
283 x_sample_windows = make_sample_windows(filter_x_scale, filter_x_support_radius, filter->func,
284 dest_width, src_width, x_scale);
285 y_sample_windows = make_sample_windows(filter_y_scale, filter_y_support_radius, filter->func,
286 dest_height, src_height, y_scale);
287
288 temp_image = (unsigned char*)malloc(num_channels * dest_width * src_height);
289
290 zoom_unidirectional(temp_image, src, num_channels, x_sample_windows,
291 dest_width, src_height,
292 num_channels, num_channels,
293 dest_row_stride, src_row_stride);
294 zoom_unidirectional(dest, temp_image, num_channels, y_sample_windows,
295 dest_height, dest_width,
296 dest_row_stride, dest_row_stride,
297 num_channels, num_channels);
298
299 free(temp_image);
300
301 free_sample_windows(x_sample_windows, dest_width);
302 free_sample_windows(y_sample_windows, dest_height);
303 }
304
305 #ifdef TEST_ZOOM
306 #include <stdio.h>
307
308 #include "readimage.h"
309 #include "writeimage.h"
310
311 int
main(int argc,char * argv[])312 main (int argc, char *argv[])
313 {
314 unsigned char *src, *dst;
315 int src_width, src_height;
316 int dst_width, dst_height;
317 void *png_write_data;
318
319 if (argc != 5)
320 {
321 fprintf(stderr, "Usage: %s <in-image> <out-width> <out-height> <out-image>\n", argv[0]);
322 return 1;
323 }
324
325 src = read_image(argv[1], &src_width, &src_height);
326 assert(src != 0);
327
328 dst_width = atoi(argv[2]);
329 dst_height = atoi(argv[3]);
330
331 dst = (unsigned char*)malloc(3 * dst_width * dst_height);
332
333 zoom_image(dst, src, get_filter(FILTER_TRIANGLE), 3,
334 dst_width, dst_height, dst_width * 3,
335 src_width, src_height, src_width * 3);
336
337 write_image(argv[4], dst_width, dst_height, dst, IMAGE_FORMAT_PNG);
338
339 return 0;
340
341 }
342 #endif
343