1 /*
2  * Image assembling routines
3  * Copyright (C) 2007-2008 Daniel Drake <dsd@gentoo.org>
4  * Copyright (C) 2013 Arseniy Lartsev <arseniy@chalmers.se>
5  * Copyright (C) 2015 Vasily Khoruzhick <anarsoul@gmail.com>
6  *
7  * This library is free software; you can redistribute it and/or
8  * modify it under the terms of the GNU Lesser General Public
9  * License as published by the Free Software Foundation; either
10  * version 2.1 of the License, or (at your option) any later version.
11  *
12  * This library is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
15  * Lesser General Public License for more details.
16  *
17  * You should have received a copy of the GNU Lesser General Public
18  * License along with this library; if not, write to the Free Software
19  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
20  */
21 
22 #define FP_COMPONENT "assembling"
23 
24 #include <errno.h>
25 #include <string.h>
26 
27 #include <libusb.h>
28 #include <glib.h>
29 
30 #include "fp_internal.h"
31 #include "assembling.h"
32 
calc_error(struct fpi_frame_asmbl_ctx * ctx,struct fpi_frame * first_frame,struct fpi_frame * second_frame,int dx,int dy)33 static unsigned int calc_error(struct fpi_frame_asmbl_ctx *ctx,
34 			       struct fpi_frame *first_frame,
35 			       struct fpi_frame *second_frame,
36 			       int dx,
37 			       int dy)
38 {
39 	unsigned int width, height;
40 	unsigned int x1, y1, x2, y2, err, i, j;
41 
42 	width = ctx->frame_width - (dx > 0 ? dx : -dx);
43 	height = ctx->frame_height - dy;
44 
45 	y1 = 0;
46 	y2 = dy;
47 	i = 0;
48 	err = 0;
49 	do {
50 		x1 = dx < 0 ? 0 : dx;
51 		x2 = dx < 0 ? -dx : 0;
52 		j = 0;
53 
54 		do {
55 			unsigned char v1, v2;
56 
57 
58 			v1 = ctx->get_pixel(ctx, first_frame, x1, y1);
59 			v2 = ctx->get_pixel(ctx, second_frame, x2, y2);
60 			err += v1 > v2 ? v1 - v2 : v2 - v1;
61 			j++;
62 			x1++;
63 			x2++;
64 
65 		} while (j < width);
66 		i++;
67 		y1++;
68 		y2++;
69 	} while (i < height);
70 
71 	/* Normalize error */
72 	err *= (ctx->frame_height * ctx->frame_width);
73 	err /= (height * width);
74 
75 	if (err == 0)
76 		return INT_MAX;
77 
78 	return err;
79 }
80 
81 /* This function is rather CPU-intensive. It's better to use hardware
82  * to detect movement direction when possible.
83  */
find_overlap(struct fpi_frame_asmbl_ctx * ctx,struct fpi_frame * first_frame,struct fpi_frame * second_frame,unsigned int * min_error)84 static void find_overlap(struct fpi_frame_asmbl_ctx *ctx,
85 			 struct fpi_frame *first_frame,
86 			 struct fpi_frame *second_frame,
87 			 unsigned int *min_error)
88 {
89 	int dx, dy;
90 	unsigned int err;
91 	*min_error = 255 * ctx->frame_height * ctx->frame_width;
92 
93 	/* Seeking in horizontal and vertical dimensions,
94 	 * for horizontal dimension we'll check only 8 pixels
95 	 * in both directions. For vertical direction diff is
96 	 * rarely less than 2, so start with it.
97 	 */
98 	for (dy = 2; dy < ctx->frame_height; dy++) {
99 		for (dx = -8; dx < 8; dx++) {
100 			err = calc_error(ctx, first_frame, second_frame,
101 				dx, dy);
102 			if (err < *min_error) {
103 				*min_error = err;
104 				second_frame->delta_x = -dx;
105 				second_frame->delta_y = dy;
106 			}
107 		}
108 	}
109 }
110 
do_movement_estimation(struct fpi_frame_asmbl_ctx * ctx,GSList * stripes,size_t num_stripes,gboolean reverse)111 static unsigned int do_movement_estimation(struct fpi_frame_asmbl_ctx *ctx,
112 			    GSList *stripes, size_t num_stripes,
113 			    gboolean reverse)
114 {
115 	GSList *list_entry = stripes;
116 	GTimer *timer;
117 	int frame = 1;
118 	struct fpi_frame *prev_stripe = list_entry->data;
119 	unsigned int min_error;
120 	/* Max error is width * height * 255, for AES2501 which has the largest
121 	 * sensor its 192*16*255 = 783360. So for 32bit value it's ~5482 frame before
122 	 * we might get int overflow. Use 64bit value here to prevent integer overflow
123 	 */
124 	unsigned long long total_error = 0;
125 
126 	list_entry = g_slist_next(list_entry);
127 
128 	timer = g_timer_new();
129 	do {
130 		struct fpi_frame *cur_stripe = list_entry->data;
131 
132 		if (reverse) {
133 			find_overlap(ctx, prev_stripe, cur_stripe, &min_error);
134 			prev_stripe->delta_y = -prev_stripe->delta_y;
135 			prev_stripe->delta_x = -prev_stripe->delta_x;
136 		}
137 		else
138 			find_overlap(ctx, cur_stripe, prev_stripe, &min_error);
139 		total_error += min_error;
140 
141 		frame++;
142 		prev_stripe = cur_stripe;
143 		list_entry = g_slist_next(list_entry);
144 
145 	} while (frame < num_stripes);
146 
147 	g_timer_stop(timer);
148 	fp_dbg("calc delta completed in %f secs", g_timer_elapsed(timer, NULL));
149 	g_timer_destroy(timer);
150 
151 	return total_error / num_stripes;
152 }
153 
fpi_do_movement_estimation(struct fpi_frame_asmbl_ctx * ctx,GSList * stripes,size_t num_stripes)154 void fpi_do_movement_estimation(struct fpi_frame_asmbl_ctx *ctx,
155 			    GSList *stripes, size_t num_stripes)
156 {
157 	int err, rev_err;
158 	err = do_movement_estimation(ctx, stripes, num_stripes, FALSE);
159 	rev_err = do_movement_estimation(ctx, stripes, num_stripes, TRUE);
160 	fp_dbg("errors: %d rev: %d", err, rev_err);
161 	if (err < rev_err) {
162 		do_movement_estimation(ctx, stripes, num_stripes, FALSE);
163 	}
164 }
165 
aes_blit_stripe(struct fpi_frame_asmbl_ctx * ctx,struct fp_img * img,struct fpi_frame * stripe,int x,int y)166 static inline void aes_blit_stripe(struct fpi_frame_asmbl_ctx *ctx,
167 				   struct fp_img *img,
168 				   struct fpi_frame *stripe,
169 				   int x, int y)
170 {
171 	unsigned int ix, iy;
172 	unsigned int fx, fy;
173 	unsigned int width, height;
174 
175 	/* Find intersection */
176 	if (x < 0) {
177 		width = ctx->frame_width + x;
178 		ix = 0;
179 		fx = -x;
180 	} else {
181 		ix = x;
182 		fx = 0;
183 		width = ctx->frame_width;
184 	}
185 	if ((ix + width) > img->width)
186 		width = img->width - ix;
187 
188 	if (y < 0) {
189 		iy = 0;
190 		fy = -y;
191 		height = ctx->frame_height + y;
192 	} else {
193 		iy = y;
194 		fy = 0;
195 		height = ctx->frame_height;
196 	}
197 
198 	if (fx > ctx->frame_width)
199 		return;
200 
201 	if (fy > ctx->frame_height)
202 		return;
203 
204 	if (ix > img->width)
205 		return;
206 
207 	if (iy > img->height)
208 		return;
209 
210 	if ((iy + height) > img->height)
211 		height = img->height - iy;
212 
213 	for (; fy < height; fy++, iy++) {
214 		if (x < 0) {
215 			ix = 0;
216 			fx = -x;
217 		} else {
218 			ix = x;
219 			fx = 0;
220 		}
221 		for (; fx < width; fx++, ix++) {
222 			img->data[ix + (iy * img->width)] = ctx->get_pixel(ctx, stripe, fx, fy);
223 		}
224 	}
225 }
226 
fpi_assemble_frames(struct fpi_frame_asmbl_ctx * ctx,GSList * stripes,size_t stripes_len)227 struct fp_img *fpi_assemble_frames(struct fpi_frame_asmbl_ctx *ctx,
228 			    GSList *stripes, size_t stripes_len)
229 {
230 	GSList *stripe;
231 	struct fp_img *img;
232 	int height = 0;
233 	int i, y, x;
234 	gboolean reverse = FALSE;
235 	struct fpi_frame *fpi_frame;
236 
237 	BUG_ON(stripes_len == 0);
238 	BUG_ON(ctx->image_width < ctx->frame_width);
239 
240 	/* Calculate height */
241 	i = 0;
242 	stripe = stripes;
243 
244 	/* No offset for 1st image */
245 	fpi_frame = stripe->data;
246 	fpi_frame->delta_x = 0;
247 	fpi_frame->delta_y = 0;
248 	do {
249 		fpi_frame = stripe->data;
250 
251 		height += fpi_frame->delta_y;
252 		i++;
253 		stripe = g_slist_next(stripe);
254 	} while (i < stripes_len);
255 
256 	fp_dbg("height is %d", height);
257 
258 	if (height < 0) {
259 		reverse = TRUE;
260 		height = -height;
261 	}
262 
263 	/* For last frame */
264 	height += ctx->frame_height;
265 
266 	/* Create buffer big enough for max image */
267 	img = fpi_img_new(ctx->image_width * height);
268 	img->flags = FP_IMG_COLORS_INVERTED;
269 	img->flags |= reverse ? 0 :  FP_IMG_H_FLIPPED | FP_IMG_V_FLIPPED;
270 	img->width = ctx->image_width;
271 	img->height = height;
272 
273 	/* Assemble stripes */
274 	i = 0;
275 	stripe = stripes;
276 	y = reverse ? (height - ctx->frame_height) : 0;
277 	x = (ctx->image_width - ctx->frame_width) / 2;
278 
279 	do {
280 		fpi_frame = stripe->data;
281 
282 		y += fpi_frame->delta_y;
283 		x += fpi_frame->delta_x;
284 
285 		aes_blit_stripe(ctx, img, fpi_frame, x, y);
286 
287 		stripe = g_slist_next(stripe);
288 		i++;
289 	} while (i < stripes_len);
290 
291 	return img;
292 }
293 
cmpint(const void * p1,const void * p2,gpointer data)294 static int cmpint(const void *p1, const void *p2, gpointer data)
295 {
296 	int a = *((int *)p1);
297 	int b = *((int *)p2);
298 	if (a < b)
299 		return -1;
300 	else if (a == b)
301 		return 0;
302 	else
303 		return 1;
304 }
305 
median_filter(int * data,int size,int filtersize)306 static void median_filter(int *data, int size, int filtersize)
307 {
308 	int i;
309 	int *result = (int *)g_malloc0(size*sizeof(int));
310 	int *sortbuf = (int *)g_malloc0(filtersize*sizeof(int));
311 	for (i = 0; i < size; i++) {
312 		int i1 = i - (filtersize-1)/2;
313 		int i2 = i + (filtersize-1)/2;
314 		if (i1 < 0)
315 			i1 = 0;
316 		if (i2 >= size)
317 			i2 = size-1;
318 		g_memmove(sortbuf, data+i1, (i2-i1+1)*sizeof(int));
319 		g_qsort_with_data(sortbuf, i2-i1+1, sizeof(int), cmpint, NULL);
320 		result[i] = sortbuf[(i2-i1+1)/2];
321 	}
322 	memmove(data, result, size*sizeof(int));
323 	g_free(result);
324 	g_free(sortbuf);
325 }
326 
interpolate_lines(struct fpi_line_asmbl_ctx * ctx,GSList * line1,float y1,GSList * line2,float y2,unsigned char * output,float yi,int size)327 static void interpolate_lines(struct fpi_line_asmbl_ctx *ctx,
328 			      GSList *line1, float y1, GSList *line2,
329 			      float y2, unsigned char *output, float yi, int size)
330 {
331 	int i;
332 	unsigned char p1, p2;
333 
334 	if (!line1 || !line2)
335 		return;
336 
337 	for (i = 0; i < size; i++) {
338 		p1 = ctx->get_pixel(ctx, line1, i);
339 		p2 = ctx->get_pixel(ctx, line2, i);
340 		output[i] = (float)p1
341 			    + (yi - y1)/(y2 - y1)*(p2 - p1);
342 	}
343 }
344 
min(int a,int b)345 static int min(int a, int b) {return (a < b) ? a : b; }
346 
347 /* Rescale image to account for variable swiping speed */
fpi_assemble_lines(struct fpi_line_asmbl_ctx * ctx,GSList * lines,size_t lines_len)348 struct fp_img *fpi_assemble_lines(struct fpi_line_asmbl_ctx *ctx,
349 				  GSList *lines, size_t lines_len)
350 {
351 	/* Number of output lines per distance between two scanners */
352 	int i;
353 	GSList *row1, *row2;
354 	float y = 0.0;
355 	int line_ind = 0;
356 	int *offsets = (int *)g_malloc0((lines_len / 2) * sizeof(int));
357 	unsigned char *output = g_malloc0(ctx->line_width * ctx->max_height);
358 	struct fp_img *img;
359 
360 	fp_dbg("%llu", g_get_real_time());
361 
362 	row1 = lines;
363 	for (i = 0; (i < lines_len - 1) && row1; i += 2) {
364 		int bestmatch = i;
365 		int bestdiff = 0;
366 		int j, firstrow, lastrow;
367 
368 		firstrow = i + 1;
369 		lastrow = min(i + ctx->max_search_offset, lines_len - 1);
370 
371 		row2 = g_slist_next(row1);
372 		for (j = firstrow; j <= lastrow; j++) {
373 			int diff = ctx->get_deviation(ctx,
374 					row1,
375 					row2);
376 			if ((j == firstrow) || (diff < bestdiff)) {
377 				bestdiff = diff;
378 				bestmatch = j;
379 			}
380 			row2 = g_slist_next(row2);
381 		}
382 		offsets[i / 2] = bestmatch - i;
383 		fp_dbg("%d", offsets[i / 2]);
384 		row1 = g_slist_next(row1);
385 		if (row1)
386 			row1 = g_slist_next(row1);
387 	}
388 
389 	median_filter(offsets, (lines_len / 2) - 1, ctx->median_filter_size);
390 
391 	fp_dbg("offsets_filtered: %llu", g_get_real_time());
392 	for (i = 0; i <= (lines_len / 2) - 1; i++)
393 		fp_dbg("%d", offsets[i]);
394 	row1 = lines;
395 	for (i = 0; i < lines_len - 1; i++, row1 = g_slist_next(row1)) {
396 		int offset = offsets[i/2];
397 		if (offset > 0) {
398 			float ynext = y + (float)ctx->resolution / offset;
399 			while (line_ind < ynext) {
400 				if (line_ind > ctx->max_height - 1)
401 					goto out;
402 				interpolate_lines(ctx,
403 					row1, y,
404 					g_slist_next(row1),
405 					ynext,
406 					output + line_ind * ctx->line_width,
407 					line_ind,
408 					ctx->line_width);
409 				line_ind++;
410 			}
411 			y = ynext;
412 		}
413 	}
414 out:
415 	img = fpi_img_new(ctx->line_width * line_ind);
416 	img->height = line_ind;
417 	img->width = ctx->line_width;
418 	img->flags = FP_IMG_V_FLIPPED;
419 	g_memmove(img->data, output, ctx->line_width * line_ind);
420 	g_free(offsets);
421 	g_free(output);
422 	return img;
423 }
424