Mercurial > hgrepos > Python2 > PyMuPDF

comparison mupdf-source/source/fitz/skew.c @ 2:b50eed0cc0ef upstream

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
ADD: MuPDF v1.26.7: the MuPDF source as downloaded by a default build of PyMuPDF 1.26.4. The directory name has changed: no version number in the expanded directory now.
author Franz Glasner <fzglas.hg@dom66.de>
date Mon, 15 Sep 2025 11:43:07 +0200
parents
children
comparison
equal deleted inserted replaced
1:1d09e1dec1d9 2:b50eed0cc0ef
1 // Copyright (C) 2004-2021 Artifex Software, Inc.
2 //
3 // This file is part of MuPDF.
4 //
5 // MuPDF is free software: you can redistribute it and/or modify it under the
6 // terms of the GNU Affero General Public License as published by the Free
7 // Software Foundation, either version 3 of the License, or (at your option)
8 // any later version.
9 //
10 // MuPDF is distributed in the hope that it will be useful, but WITHOUT ANY
11 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
12 // FOR A PARTICULAR PURPOSE. See the GNU Affero General Public License for more
13 // details.
14 //
15 // You should have received a copy of the GNU Affero General Public License
16 // along with MuPDF. If not, see <https://www.gnu.org/licenses/agpl-3.0.en.html>
17 //
18 // Alternative licensing terms are available from the licensor.
19 // For commercial licensing, see <https://www.artifex.com/> or contact
20 // Artifex Software, Inc., 39 Mesa Street, Suite 108A, San Francisco,
21 // CA 94129, USA, for further information.
22
23 #include "mupdf/fitz.h"
24
25 #include <math.h>
26 #include <assert.h>
27 #include <limits.h>
28
29 #if ARCH_HAS_SSE
30 #include <emmintrin.h>
31 #include <smmintrin.h>
32 #endif
33
34 /* Uncomment the following to enable Debugging printfs. */
35 /* #define DEBUG_PRINT_WORKING */
36
37 enum {
38 NUM_SKEW_COLS = 4,
39 SKEW_COL_OFFSET = 96,
40 SKEW_COL_WIDTH = 96,
41 SKEW_MAX_DIFF = 16
42 };
43
44 typedef struct
45 {
46 void *src;
47 uint32_t src_w;
48 uint32_t src_h;
49 uint32_t y;
50 uint32_t offsets[NUM_SKEW_COLS*2];
51 uint16_t *tables[NUM_SKEW_COLS*2];
52 uint32_t corr_height;
53 } fz_skew;
54
55 /* Finalise a rescaler instance. */
56 static void
57 fz_drop_skew(fz_context *ctx, fz_skew *skew)
58 {
59 if (skew == NULL)
60 return;
61
62 fz_free(ctx, skew->tables[0]);
63 fz_free(ctx, skew);
64 }
65
66 /* Initialise a skew instance. */
67 static fz_skew *
68 fz_new_skew(fz_context *ctx,
69 unsigned int src_w,
70 unsigned int src_h)
71 {
72 fz_skew *skew = fz_malloc_struct(ctx, fz_skew);
73 int i;
74
75 skew->src_w = src_w;
76 skew->src_h = src_h;
77 skew->y = 0;
78 skew->corr_height = src_h-100; /* FIXME */
79 fz_try(ctx)
80 skew->tables[0] = (uint16_t *)fz_malloc(ctx, skew->src_h * sizeof(uint16_t) * NUM_SKEW_COLS * 2);
81 fz_catch(ctx)
82 {
83 fz_drop_skew(ctx, skew);
84 fz_rethrow(ctx);
85 }
86
87 for (i = 1; i < NUM_SKEW_COLS * 2; i++)
88 skew->tables[i] = skew->tables[0] + i*skew->src_h;
89
90 for (i = 0; i < NUM_SKEW_COLS; i++)
91 {
92 int offset = skew->src_w * (i+1) / (NUM_SKEW_COLS+1) - SKEW_COL_WIDTH/2;
93 skew->offsets[2*i ] = offset - SKEW_COL_OFFSET;
94 skew->offsets[2*i+1] = offset + SKEW_COL_OFFSET;
95 }
96
97 return skew;
98 }
99
100 static int sum_c(const uint8_t *data)
101 {
102 int i;
103 uint32_t sum = 0;
104 for (i = SKEW_COL_WIDTH; i > 0; i--)
105 sum += *data++;
106 return sum;
107 }
108
109 #if ARCH_HAS_SSE
110 static inline int sum_sse(const uint8_t *data)
111 {
112 __m128i mm0, mm1, mm2, mm3, mm4, mm5;
113 __m128i mm_zero = _mm_set1_epi32(0);
114
115 mm0 = _mm_loadu_si128((const __m128i *)(data )); // mm0 = ppoonnmmllkkjjiihhggffeeddccbbaa
116 mm1 = _mm_loadu_si128((const __m128i *)(data+16)); // mm1 = ppoonnmmllkkjjiihhggffeeddccbbaa
117 mm2 = _mm_loadu_si128((const __m128i *)(data+32)); // mm2 = ppoonnmmllkkjjiihhggffeeddccbbaa
118 mm3 = _mm_loadu_si128((const __m128i *)(data+48)); // mm3 = ppoonnmmllkkjjiihhggffeeddccbbaa
119 mm4 = _mm_loadu_si128((const __m128i *)(data+64)); // mm4 = ppoonnmmllkkjjiihhggffeeddccbbaa
120 mm5 = _mm_loadu_si128((const __m128i *)(data+80)); // mm5 = ppoonnmmllkkjjiihhggffeeddccbbaa
121
122 mm0 = _mm_sad_epu8(mm0, mm_zero); // Max value in each half is 8*255
123 mm1 = _mm_sad_epu8(mm1, mm_zero);
124 mm2 = _mm_sad_epu8(mm2, mm_zero);
125 mm3 = _mm_sad_epu8(mm3, mm_zero);
126 mm4 = _mm_sad_epu8(mm4, mm_zero);
127 mm5 = _mm_sad_epu8(mm5, mm_zero);
128
129 mm0 = _mm_add_epi64(mm0, mm1); // Max value in each half is 2*8*255
130 mm2 = _mm_add_epi64(mm2, mm3); // Max value in each half is 2*8*255
131 mm4 = _mm_add_epi64(mm4, mm5); // Max value in each half is 2*8*255
132 mm0 = _mm_add_epi64(mm0, mm2); // Max value in each half is 4*8*255
133 mm0 = _mm_add_epi64(mm0, mm4); // Max value in each half is 6*8*255
134
135 mm0 = _mm_shuffle_epi32(mm0, (2<<2)+0); // Max value in each half is 10*8*255
136 mm0 = _mm_hadd_epi32(mm0, mm0); // Max value in bottom bits is 20*8*255 - still fits in an unsigned 16bit word.
137
138 return _mm_extract_epi16(mm0, 0);
139 }
140 #endif
141
142 /* Process data: */
143 static void
144 fz_skew_process(fz_context *ctx, fz_skew *skew, const uint8_t *input)
145 {
146 int i;
147 int off = skew->y++;
148
149 #if ARCH_HAS_SSE
150 for (i = 0; i < NUM_SKEW_COLS * 2; i++)
151 skew->tables[i][off] = sum_sse(input + skew->offsets[i]);
152 #else
153 for (i = 0; i < NUM_SKEW_COLS * 2; i++)
154 skew->tables[i][off] = sum_c(input + skew->offsets[i]);
155 #endif
156
157 /* Some debug code; if enabled this writes the summed value back
158 * in to give us a visual indication of where we are looking for
159 * correspondences. */
160 #if 0
161 for (i = 0; i < NUM_SKEW_COLS * 2; i++) {
162 int v = (skew->tables[i][off] + (SKEW_COL_WIDTH/2) ) / SKEW_COL_WIDTH;
163 memset(input + skew->offsets[i], v, SKEW_COL_WIDTH);
164 }
165 #endif
166 }
167
168 static double
169 do_detect_skew(fz_context *ctx, fz_skew *skew)
170 {
171 int i, j, h, o, max_at;
172 int64_t max_sum, corr_at, corr_sum, avg_sum;
173 float ang;
174 int64_t avg = SKEW_COL_WIDTH * 255/2;
175
176 if (skew == NULL)
177 return 0;
178
179 h = skew->corr_height - 2 * SKEW_MAX_DIFF;
180
181 corr_at = 0;
182 corr_sum = 0;
183 for (i = 0; i < NUM_SKEW_COLS; i++)
184 {
185 max_at = 9999;
186 max_sum = 0;
187 avg_sum = 0;
188 for (o = -SKEW_MAX_DIFF; o <= SKEW_MAX_DIFF; o++)
189 {
190 uint16_t *t0 = skew->tables[2*i] + SKEW_MAX_DIFF;
191 uint16_t *t1 = skew->tables[2*i+1] + SKEW_MAX_DIFF + o;
192 int64_t sum = 0;
193 for (j = 0; j < h; j++)
194 sum += ((int64_t)t0[j]-avg) * ((int64_t)t1[j]-avg);
195 if (max_sum < sum)
196 max_sum = sum, max_at = o;
197 avg_sum += sum;
198 #ifdef DEBUG_PRINT_WORKING
199 printf("col %d, offset %d -> %llx\n", i, o, sum);
200 #endif
201 }
202 avg_sum /= (SKEW_MAX_DIFF+1)*2;
203 #ifdef DEBUG_PRINT_WORKING
204 ang = (180.0 / 3.1415942) * atan(max_at / (double)(SKEW_COL_OFFSET * 2));
205 printf("col %d max: offset %d -> %llx ang=%g\n", i, max_at, max_sum, ang);
206 #endif
207 /* Subtract the average from the maximum; we judge the significance of a
208 * match by how far it exceeds the average. max_sum becomes 'significance'. */
209 max_sum -= avg_sum;
210 #ifdef DEBUG_PRINT_WORKING
211 printf("Significance: %llx\n", max_sum - avg_sum);
212 #endif
213 corr_at += max_at * max_sum;
214 corr_sum += max_sum;
215 }
216 ang = (180.0 / 3.1415942) * atan(corr_at / (double)(corr_sum * SKEW_COL_OFFSET * 2));
217
218 if (ang < -45 || ang > 45)
219 return 0;
220
221 return ang;
222 }
223
224 double fz_detect_skew(fz_context *ctx, fz_pixmap *pix)
225 {
226 fz_skew *skew = fz_new_skew(ctx, pix->w, pix->h);
227 int y;
228 uint8_t *ptr;
229 ptrdiff_t stride;
230 double angle;
231 fz_pixmap *pix2 = NULL;
232
233 fz_var(pix2);
234
235 fz_try(ctx)
236 {
237 if (pix->n != 1)
238 {
239 pix2 = fz_convert_pixmap(ctx, pix, fz_device_gray(ctx), NULL, NULL, fz_default_color_params, 0);
240 ptr = pix2->samples;
241 stride = pix2->stride;
242 }
243 else
244 {
245 ptr = pix->samples;
246 stride = pix->stride;
247 }
248
249 for (y = pix->h; y > 0; y--)
250 {
251 fz_skew_process(ctx, skew, ptr);
252 ptr += stride;
253 }
254
255 angle = do_detect_skew(ctx, skew);
256 }
257 fz_always(ctx)
258 {
259 fz_drop_pixmap(ctx, pix2);
260 fz_drop_skew(ctx, skew);
261 }
262 fz_catch(ctx)
263 fz_rethrow(ctx);
264
265 return angle;
266 }