Coverage Report

Created: 2025-06-13 07:15

/src/leptonica/src/projective.c
Line
Count
Source (jump to first uncovered line)
1
/*====================================================================*
2
 -  Copyright (C) 2001 Leptonica.  All rights reserved.
3
 -
4
 -  Redistribution and use in source and binary forms, with or without
5
 -  modification, are permitted provided that the following conditions
6
 -  are met:
7
 -  1. Redistributions of source code must retain the above copyright
8
 -     notice, this list of conditions and the following disclaimer.
9
 -  2. Redistributions in binary form must reproduce the above
10
 -     copyright notice, this list of conditions and the following
11
 -     disclaimer in the documentation and/or other materials
12
 -     provided with the distribution.
13
 -
14
 -  THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
15
 -  ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
16
 -  LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
17
 -  A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL ANY
18
 -  CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
19
 -  EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
20
 -  PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
21
 -  PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY
22
 -  OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
23
 -  NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
24
 -  SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
25
 *====================================================================*/
26
27
/*!
28
 * \file projective.c
29
 * <pre>
30
 *
31
 *      Projective (4 pt) image transformation using a sampled
32
 *      (to nearest integer) transform on each dest point
33
 *           PIX      *pixProjectiveSampledPta()
34
 *           PIX      *pixProjectiveSampled()
35
 *
36
 *      Projective (4 pt) image transformation using interpolation
37
 *      (or area mapping) for anti-aliasing images that are
38
 *      2, 4, or 8 bpp gray, or colormapped, or 32 bpp RGB
39
 *           PIX      *pixProjectivePta()
40
 *           PIX      *pixProjective()
41
 *           PIX      *pixProjectivePtaColor()
42
 *           PIX      *pixProjectiveColor()
43
 *           PIX      *pixProjectivePtaGray()
44
 *           PIX      *pixProjectiveGray()
45
 *
46
 *      Projective transform including alpha (blend) component
47
 *           PIX      *pixProjectivePtaWithAlpha()
48
 *
49
 *      Projective coordinate transformation
50
 *           l_int32   getProjectiveXformCoeffs()
51
 *           l_int32   projectiveXformSampledPt()
52
 *           l_int32   projectiveXformPt()
53
 *
54
 *      A projective transform can be specified as a specific functional
55
 *      mapping between 4 points in the source and 4 points in the dest.
56
 *      It preserves straight lines, but is less stable than a bilinear
57
 *      transform, because it contains a division by a quantity that
58
 *      can get arbitrarily small.)
59
 *
60
 *      We give both a projective coordinate transformation and
61
 *      two projective image transformations.
62
 *
63
 *      For the former, we ask for the coordinate value (x',y')
64
 *      in the transformed space for any point (x,y) in the original
65
 *      space.  The coefficients of the transformation are found by
66
 *      solving 8 simultaneous equations for the 8 coordinates of
67
 *      the 4 points in src and dest.  The transformation can then
68
 *      be used to compute the associated image transform, by
69
 *      computing, for each dest pixel, the relevant pixel(s) in
70
 *      the source.  This can be done either by taking the closest
71
 *      src pixel to each transformed dest pixel ("sampling") or
72
 *      by doing an interpolation and averaging over 4 source
73
 *      pixels with appropriate weightings ("interpolated").
74
 *
75
 *      A typical application would be to remove keystoning
76
 *      due to a projective transform in the imaging system.
77
 *
78
 *      The projective transform is given by specifying two equations:
79
 *
80
 *          x' = (ax + by + c) / (gx + hy + 1)
81
 *          y' = (dx + ey + f) / (gx + hy + 1)
82
 *
83
 *      where the eight coefficients have been computed from four
84
 *      sets of these equations, each for two corresponding data pts.
85
 *      In practice, once the coefficients are known, we use the
86
 *      equations "backwards": for each point (x,y) in the dest image,
87
 *      these two equations are used to compute the corresponding point
88
 *      (x',y') in the src.  That computed point in the src is then used
89
 *      to determine the corresponding dest pixel value in one of two ways:
90
 *
91
 *       ~ sampling: simply take the value of the src pixel in which this
92
 *                   point falls
93
 *       ~ interpolation: take appropriate linear combinations of the
94
 *                        four src pixels that this dest pixel would
95
 *                        overlap, with the coefficients proportional
96
 *                        to the amount of overlap
97
 *
98
 *      For small warp where there is little scale change, (e.g.,
99
 *      for rotation) area mapping is nearly equivalent to interpolation.
100
 *
101
 *      Typical relative timing of pointwise transforms (sampled = 1.0):
102
 *      8 bpp:   sampled        1.0
103
 *               interpolated   1.5
104
 *      32 bpp:  sampled        1.0
105
 *               interpolated   1.6
106
 *      Additionally, the computation time/pixel is nearly the same
107
 *      for 8 bpp and 32 bpp, for both sampled and interpolated.
108
 * </pre>
109
 */
110
111
#ifdef HAVE_CONFIG_H
112
#include <config_auto.h>
113
#endif  /* HAVE_CONFIG_H */
114
115
#include <string.h>
116
#include <math.h>
117
#include "allheaders.h"
118
119
extern l_float32  AlphaMaskBorderVals[2];
120
121
/*------------------------------------------------------------n
122
 *            Sampled projective image transformation          *
123
 *-------------------------------------------------------------*/
124
/*!
125
 * \brief   pixProjectiveSampledPta()
126
 *
127
 * \param[in]    pixs      all depths
128
 * \param[in]    ptad      4 pts of final coordinate space
129
 * \param[in]    ptas      4 pts of initial coordinate space
130
 * \param[in]    incolor   L_BRING_IN_WHITE, L_BRING_IN_BLACK
131
 * \return  pixd, or NULL on error
132
 *
133
 * <pre>
134
 * Notes:
135
 *      (1) Brings in either black or white pixels from the boundary.
136
 *      (2) Retains colormap, which you can do for a sampled transform..
137
 *      (3) No 3 of the 4 points may be collinear.
138
 *      (4) For 8 and 32 bpp pix, better quality is obtained by the
139
 *          somewhat slower pixProjectivePta().  See that
140
 *          function for relative timings between sampled and interpolated.
141
 * </pre>
142
 */
143
PIX *
144
pixProjectiveSampledPta(PIX     *pixs,
145
                        PTA     *ptad,
146
                        PTA     *ptas,
147
                        l_int32  incolor)
148
0
{
149
0
l_float32  *vc;
150
0
PIX        *pixd;
151
152
0
    if (!pixs)
153
0
        return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
154
0
    if (!ptas)
155
0
        return (PIX *)ERROR_PTR("ptas not defined", __func__, NULL);
156
0
    if (!ptad)
157
0
        return (PIX *)ERROR_PTR("ptad not defined", __func__, NULL);
158
0
    if (incolor != L_BRING_IN_WHITE && incolor != L_BRING_IN_BLACK)
159
0
        return (PIX *)ERROR_PTR("invalid incolor", __func__, NULL);
160
0
    if (ptaGetCount(ptas) != 4)
161
0
        return (PIX *)ERROR_PTR("ptas count not 4", __func__, NULL);
162
0
    if (ptaGetCount(ptad) != 4)
163
0
        return (PIX *)ERROR_PTR("ptad count not 4", __func__, NULL);
164
165
        /* Get backwards transform from dest to src, and apply it */
166
0
    getProjectiveXformCoeffs(ptad, ptas, &vc);
167
0
    pixd = pixProjectiveSampled(pixs, vc, incolor);
168
0
    LEPT_FREE(vc);
169
170
0
    return pixd;
171
0
}
172
173
174
/*!
175
 * \brief   pixProjectiveSampled()
176
 *
177
 * \param[in]    pixs      all depths
178
 * \param[in]    vc        vector of 8 coefficients for projective transform
179
 * \param[in]    incolor   L_BRING_IN_WHITE, L_BRING_IN_BLACK
180
 * \return  pixd, or NULL on error
181
 *
182
 * <pre>
183
 * Notes:
184
 *      (1) Brings in either black or white pixels from the boundary.
185
 *      (2) Retains colormap, which you can do for a sampled transform..
186
 *      (3) For 8 or 32 bpp, much better quality is obtained by the
187
 *          somewhat slower pixProjective().  See that function
188
 *          for relative timings between sampled and interpolated.
189
 * </pre>
190
 */
191
PIX *
192
pixProjectiveSampled(PIX        *pixs,
193
                     l_float32  *vc,
194
                     l_int32     incolor)
195
0
{
196
0
l_int32     i, j, w, h, d, x, y, wpls, wpld, color, cmapindex;
197
0
l_uint32    val;
198
0
l_uint32   *datas, *datad, *lines, *lined;
199
0
PIX        *pixd;
200
0
PIXCMAP    *cmap;
201
202
0
    if (!pixs)
203
0
        return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
204
0
    if (!vc)
205
0
        return (PIX *)ERROR_PTR("vc not defined", __func__, NULL);
206
0
    if (incolor != L_BRING_IN_WHITE && incolor != L_BRING_IN_BLACK)
207
0
        return (PIX *)ERROR_PTR("invalid incolor", __func__, NULL);
208
0
    pixGetDimensions(pixs, &w, &h, &d);
209
0
    if (d != 1 && d != 2 && d != 4 && d != 8 && d != 32)
210
0
        return (PIX *)ERROR_PTR("depth not 1, 2, 4, 8 or 16", __func__, NULL);
211
212
        /* Init all dest pixels to color to be brought in from outside */
213
0
    pixd = pixCreateTemplate(pixs);
214
0
    if ((cmap = pixGetColormap(pixs)) != NULL) {
215
0
        if (incolor == L_BRING_IN_WHITE)
216
0
            color = 1;
217
0
        else
218
0
            color = 0;
219
0
        pixcmapAddBlackOrWhite(cmap, color, &cmapindex);
220
0
        pixSetAllArbitrary(pixd, cmapindex);
221
0
    } else {
222
0
        if ((d == 1 && incolor == L_BRING_IN_WHITE) ||
223
0
            (d > 1 && incolor == L_BRING_IN_BLACK)) {
224
0
            pixClearAll(pixd);
225
0
        } else {
226
0
            pixSetAll(pixd);
227
0
        }
228
0
    }
229
230
        /* Scan over the dest pixels */
231
0
    datas = pixGetData(pixs);
232
0
    wpls = pixGetWpl(pixs);
233
0
    datad = pixGetData(pixd);
234
0
    wpld = pixGetWpl(pixd);
235
0
    for (i = 0; i < h; i++) {
236
0
        lined = datad + i * wpld;
237
0
        for (j = 0; j < w; j++) {
238
0
            projectiveXformSampledPt(vc, j, i, &x, &y);
239
0
            if (x < 0 || y < 0 || x >=w || y >= h)
240
0
                continue;
241
0
            lines = datas + y * wpls;
242
0
            if (d == 1) {
243
0
                val = GET_DATA_BIT(lines, x);
244
0
                SET_DATA_BIT_VAL(lined, j, val);
245
0
            } else if (d == 8) {
246
0
                val = GET_DATA_BYTE(lines, x);
247
0
                SET_DATA_BYTE(lined, j, val);
248
0
            } else if (d == 32) {
249
0
                lined[j] = lines[x];
250
0
            } else if (d == 2) {
251
0
                val = GET_DATA_DIBIT(lines, x);
252
0
                SET_DATA_DIBIT(lined, j, val);
253
0
            } else if (d == 4) {
254
0
                val = GET_DATA_QBIT(lines, x);
255
0
                SET_DATA_QBIT(lined, j, val);
256
0
            }
257
0
        }
258
0
    }
259
260
0
    return pixd;
261
0
}
262
263
264
/*---------------------------------------------------------------------*
265
 *            Interpolated projective image transformation             *
266
 *---------------------------------------------------------------------*/
267
/*!
268
 * \brief   pixProjectivePta()
269
 *
270
 * \param[in]    pixs      all depths; colormap ok
271
 * \param[in]    ptad      4 pts of final coordinate space
272
 * \param[in]    ptas      4 pts of initial coordinate space
273
 * \param[in]    incolor   L_BRING_IN_WHITE, L_BRING_IN_BLACK
274
 * \return  pixd, or NULL on error
275
 *
276
 * <pre>
277
 * Notes:
278
 *      (1) Brings in either black or white pixels from the boundary
279
 *      (2) Removes any existing colormap, if necessary, before transforming
280
 * </pre>
281
 */
282
PIX *
283
pixProjectivePta(PIX     *pixs,
284
                 PTA     *ptad,
285
                 PTA     *ptas,
286
                 l_int32  incolor)
287
0
{
288
0
l_int32   d;
289
0
l_uint32  colorval;
290
0
PIX      *pixt1, *pixt2, *pixd;
291
292
0
    if (!pixs)
293
0
        return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
294
0
    if (!ptas)
295
0
        return (PIX *)ERROR_PTR("ptas not defined", __func__, NULL);
296
0
    if (!ptad)
297
0
        return (PIX *)ERROR_PTR("ptad not defined", __func__, NULL);
298
0
    if (incolor != L_BRING_IN_WHITE && incolor != L_BRING_IN_BLACK)
299
0
        return (PIX *)ERROR_PTR("invalid incolor", __func__, NULL);
300
0
    if (ptaGetCount(ptas) != 4)
301
0
        return (PIX *)ERROR_PTR("ptas count not 4", __func__, NULL);
302
0
    if (ptaGetCount(ptad) != 4)
303
0
        return (PIX *)ERROR_PTR("ptad count not 4", __func__, NULL);
304
305
0
    if (pixGetDepth(pixs) == 1)
306
0
        return pixProjectiveSampledPta(pixs, ptad, ptas, incolor);
307
308
        /* Remove cmap if it exists, and unpack to 8 bpp if necessary */
309
0
    pixt1 = pixRemoveColormap(pixs, REMOVE_CMAP_BASED_ON_SRC);
310
0
    d = pixGetDepth(pixt1);
311
0
    if (d < 8)
312
0
        pixt2 = pixConvertTo8(pixt1, FALSE);
313
0
    else
314
0
        pixt2 = pixClone(pixt1);
315
0
    d = pixGetDepth(pixt2);
316
317
        /* Compute actual color to bring in from edges */
318
0
    colorval = 0;
319
0
    if (incolor == L_BRING_IN_WHITE) {
320
0
        if (d == 8)
321
0
            colorval = 255;
322
0
        else  /* d == 32 */
323
0
            colorval = 0xffffff00;
324
0
    }
325
326
0
    if (d == 8)
327
0
        pixd = pixProjectivePtaGray(pixt2, ptad, ptas, colorval);
328
0
    else  /* d == 32 */
329
0
        pixd = pixProjectivePtaColor(pixt2, ptad, ptas, colorval);
330
0
    pixDestroy(&pixt1);
331
0
    pixDestroy(&pixt2);
332
0
    return pixd;
333
0
}
334
335
336
/*!
337
 * \brief   pixProjective()
338
 *
339
 * \param[in]    pixs      all depths; colormap ok
340
 * \param[in]    vc        vector of 8 coefficients for projective transform
341
 * \param[in]    incolor   L_BRING_IN_WHITE, L_BRING_IN_BLACK
342
 * \return  pixd, or NULL on error
343
 *
344
 * <pre>
345
 * Notes:
346
 *      (1) Brings in either black or white pixels from the boundary
347
 *      (2) Removes any existing colormap, if necessary, before transforming
348
 * </pre>
349
 */
350
PIX *
351
pixProjective(PIX        *pixs,
352
              l_float32  *vc,
353
              l_int32     incolor)
354
0
{
355
0
l_int32   d;
356
0
l_uint32  colorval;
357
0
PIX      *pixt1, *pixt2, *pixd;
358
359
0
    if (!pixs)
360
0
        return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
361
0
    if (!vc)
362
0
        return (PIX *)ERROR_PTR("vc not defined", __func__, NULL);
363
364
0
    if (pixGetDepth(pixs) == 1)
365
0
        return pixProjectiveSampled(pixs, vc, incolor);
366
367
        /* Remove cmap if it exists, and unpack to 8 bpp if necessary */
368
0
    pixt1 = pixRemoveColormap(pixs, REMOVE_CMAP_BASED_ON_SRC);
369
0
    d = pixGetDepth(pixt1);
370
0
    if (d < 8)
371
0
        pixt2 = pixConvertTo8(pixt1, FALSE);
372
0
    else
373
0
        pixt2 = pixClone(pixt1);
374
0
    d = pixGetDepth(pixt2);
375
376
        /* Compute actual color to bring in from edges */
377
0
    colorval = 0;
378
0
    if (incolor == L_BRING_IN_WHITE) {
379
0
        if (d == 8)
380
0
            colorval = 255;
381
0
        else  /* d == 32 */
382
0
            colorval = 0xffffff00;
383
0
    }
384
385
0
    if (d == 8)
386
0
        pixd = pixProjectiveGray(pixt2, vc, colorval);
387
0
    else  /* d == 32 */
388
0
        pixd = pixProjectiveColor(pixt2, vc, colorval);
389
0
    pixDestroy(&pixt1);
390
0
    pixDestroy(&pixt2);
391
0
    return pixd;
392
0
}
393
394
395
/*!
396
 * \brief   pixProjectivePtaColor()
397
 *
398
 * \param[in]    pixs 32 bpp
399
 * \param[in]    ptad  4 pts of final coordinate space
400
 * \param[in]    ptas  4 pts of initial coordinate space
401
 * \param[in]    colorval e.g., 0 to bring in BLACK, 0xffffff00 for WHITE
402
 * \return  pixd, or NULL on error
403
 */
404
PIX *
405
pixProjectivePtaColor(PIX      *pixs,
406
                      PTA      *ptad,
407
                      PTA      *ptas,
408
                      l_uint32  colorval)
409
0
{
410
0
l_float32  *vc;
411
0
PIX        *pixd;
412
413
0
    if (!pixs)
414
0
        return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
415
0
    if (!ptas)
416
0
        return (PIX *)ERROR_PTR("ptas not defined", __func__, NULL);
417
0
    if (!ptad)
418
0
        return (PIX *)ERROR_PTR("ptad not defined", __func__, NULL);
419
0
    if (pixGetDepth(pixs) != 32)
420
0
        return (PIX *)ERROR_PTR("pixs must be 32 bpp", __func__, NULL);
421
0
    if (ptaGetCount(ptas) != 4)
422
0
        return (PIX *)ERROR_PTR("ptas count not 4", __func__, NULL);
423
0
    if (ptaGetCount(ptad) != 4)
424
0
        return (PIX *)ERROR_PTR("ptad count not 4", __func__, NULL);
425
426
        /* Get backwards transform from dest to src, and apply it */
427
0
    getProjectiveXformCoeffs(ptad, ptas, &vc);
428
0
    pixd = pixProjectiveColor(pixs, vc, colorval);
429
0
    LEPT_FREE(vc);
430
431
0
    return pixd;
432
0
}
433
434
435
/*!
436
 * \brief   pixProjectiveColor()
437
 *
438
 * \param[in]    pixs       32 bpp
439
 * \param[in]    vc         vector of 8 coefficients for projective transform
440
 * \param[in]    colorval   e.g., 0 to bring in BLACK, 0xffffff00 for WHITE
441
 * \return  pixd, or NULL on error
442
 */
443
PIX *
444
pixProjectiveColor(PIX        *pixs,
445
                   l_float32  *vc,
446
                   l_uint32    colorval)
447
0
{
448
0
l_int32    i, j, w, h, d, wpls, wpld;
449
0
l_uint32   val;
450
0
l_uint32  *datas, *datad, *lined;
451
0
l_float32  x, y;
452
0
PIX       *pix1, *pix2, *pixd;
453
454
0
    if (!pixs)
455
0
        return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
456
0
    pixGetDimensions(pixs, &w, &h, &d);
457
0
    if (d != 32)
458
0
        return (PIX *)ERROR_PTR("pixs must be 32 bpp", __func__, NULL);
459
0
    if (!vc)
460
0
        return (PIX *)ERROR_PTR("vc not defined", __func__, NULL);
461
462
0
    datas = pixGetData(pixs);
463
0
    wpls = pixGetWpl(pixs);
464
0
    pixd = pixCreateTemplate(pixs);
465
0
    pixSetAllArbitrary(pixd, colorval);
466
0
    datad = pixGetData(pixd);
467
0
    wpld = pixGetWpl(pixd);
468
469
        /* Iterate over destination pixels */
470
0
    for (i = 0; i < h; i++) {
471
0
        lined = datad + i * wpld;
472
0
        for (j = 0; j < w; j++) {
473
                /* Compute float src pixel location corresponding to (i,j) */
474
0
            projectiveXformPt(vc, j, i, &x, &y);
475
0
            linearInterpolatePixelColor(datas, wpls, w, h, x, y, colorval,
476
0
                                        &val);
477
0
            *(lined + j) = val;
478
0
        }
479
0
    }
480
481
        /* If rgba, transform the pixs alpha channel and insert in pixd */
482
0
    if (pixGetSpp(pixs) == 4) {
483
0
        pix1 = pixGetRGBComponent(pixs, L_ALPHA_CHANNEL);
484
0
        pix2 = pixProjectiveGray(pix1, vc, 255);  /* bring in opaque */
485
0
        pixSetRGBComponent(pixd, pix2, L_ALPHA_CHANNEL);
486
0
        pixDestroy(&pix1);
487
0
        pixDestroy(&pix2);
488
0
    }
489
490
0
    return pixd;
491
0
}
492
493
494
/*!
495
 * \brief   pixProjectivePtaGray()
496
 *
497
 * \param[in]    pixs      8 bpp
498
 * \param[in]    ptad      4 pts of final coordinate space
499
 * \param[in]    ptas      4 pts of initial coordinate space
500
 * \param[in]    grayval   0 to bring in BLACK, 255 for WHITE
501
 * \return  pixd, or NULL on error
502
 */
503
PIX *
504
pixProjectivePtaGray(PIX     *pixs,
505
                     PTA     *ptad,
506
                     PTA     *ptas,
507
                     l_uint8  grayval)
508
0
{
509
0
l_float32  *vc;
510
0
PIX        *pixd;
511
512
0
    if (!pixs)
513
0
        return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
514
0
    if (!ptas)
515
0
        return (PIX *)ERROR_PTR("ptas not defined", __func__, NULL);
516
0
    if (!ptad)
517
0
        return (PIX *)ERROR_PTR("ptad not defined", __func__, NULL);
518
0
    if (pixGetDepth(pixs) != 8)
519
0
        return (PIX *)ERROR_PTR("pixs must be 8 bpp", __func__, NULL);
520
0
    if (ptaGetCount(ptas) != 4)
521
0
        return (PIX *)ERROR_PTR("ptas count not 4", __func__, NULL);
522
0
    if (ptaGetCount(ptad) != 4)
523
0
        return (PIX *)ERROR_PTR("ptad count not 4", __func__, NULL);
524
525
        /* Get backwards transform from dest to src, and apply it */
526
0
    getProjectiveXformCoeffs(ptad, ptas, &vc);
527
0
    pixd = pixProjectiveGray(pixs, vc, grayval);
528
0
    LEPT_FREE(vc);
529
530
0
    return pixd;
531
0
}
532
533
534
535
/*!
536
 * \brief   pixProjectiveGray()
537
 *
538
 * \param[in]    pixs      8 bpp
539
 * \param[in]    vc        vector of 8 coefficients for projective transform
540
 * \param[in]    grayval   0 to bring in BLACK, 255 for WHITE
541
 * \return  pixd, or NULL on error
542
 */
543
PIX *
544
pixProjectiveGray(PIX        *pixs,
545
                  l_float32  *vc,
546
                  l_uint8     grayval)
547
0
{
548
0
l_int32    i, j, w, h, wpls, wpld, val;
549
0
l_uint32  *datas, *datad, *lined;
550
0
l_float32  x, y;
551
0
PIX       *pixd;
552
553
0
    if (!pixs)
554
0
        return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
555
0
    pixGetDimensions(pixs, &w, &h, NULL);
556
0
    if (pixGetDepth(pixs) != 8)
557
0
        return (PIX *)ERROR_PTR("pixs must be 8 bpp", __func__, NULL);
558
0
    if (!vc)
559
0
        return (PIX *)ERROR_PTR("vc not defined", __func__, NULL);
560
561
0
    datas = pixGetData(pixs);
562
0
    wpls = pixGetWpl(pixs);
563
0
    pixd = pixCreateTemplate(pixs);
564
0
    pixSetAllArbitrary(pixd, grayval);
565
0
    datad = pixGetData(pixd);
566
0
    wpld = pixGetWpl(pixd);
567
568
        /* Iterate over destination pixels */
569
0
    for (i = 0; i < h; i++) {
570
0
        lined = datad + i * wpld;
571
0
        for (j = 0; j < w; j++) {
572
                /* Compute float src pixel location corresponding to (i,j) */
573
0
            projectiveXformPt(vc, j, i, &x, &y);
574
0
            linearInterpolatePixelGray(datas, wpls, w, h, x, y, grayval, &val);
575
0
            SET_DATA_BYTE(lined, j, val);
576
0
        }
577
0
    }
578
579
0
    return pixd;
580
0
}
581
582
583
/*---------------------------------------------------------------------------*
584
 *            Projective transform including alpha (blend) component         *
585
 *---------------------------------------------------------------------------*/
586
/*!
587
 * \brief   pixProjectivePtaWithAlpha()
588
 *
589
 * \param[in]    pixs     32 bpp rgb
590
 * \param[in]    ptad     4 pts of final coordinate space
591
 * \param[in]    ptas     4 pts of initial coordinate space
592
 * \param[in]    pixg     [optional] 8 bpp, for alpha channel, can be null
593
 * \param[in]    fract    between 0.0 and 1.0, with 0.0 fully transparent
594
 *                        and 1.0 fully opaque
595
 * \param[in]    border   of pixels added to capture transformed source pixels
596
 * \return  pixd, or NULL on error
597
 *
598
 * <pre>
599
 * Notes:
600
 *      (1) The alpha channel is transformed separately from pixs,
601
 *          and aligns with it, being fully transparent outside the
602
 *          boundary of the transformed pixs.  For pixels that are fully
603
 *          transparent, a blending function like pixBlendWithGrayMask()
604
 *          will give zero weight to corresponding pixels in pixs.
605
 *      (2) If pixg is NULL, it is generated as an alpha layer that is
606
 *          partially opaque, using %fract.  Otherwise, it is cropped
607
 *          to pixs if required and %fract is ignored.  The alpha channel
608
 *          in pixs is never used.
609
 *      (3) Colormaps are removed.
610
 *      (4) When pixs is transformed, it doesn't matter what color is brought
611
 *          in because the alpha channel will be transparent (0) there.
612
 *      (5) To avoid losing source pixels in the destination, it may be
613
 *          necessary to add a border to the source pix before doing
614
 *          the projective transformation.  This can be any non-negative
615
 *          number.
616
 *      (6) The input %ptad and %ptas are in a coordinate space before
617
 *          the border is added.  Internally, we compensate for this
618
 *          before doing the projective transform on the image after
619
 *          the border is added.
620
 *      (7) The default setting for the border values in the alpha channel
621
 *          is 0 (transparent) for the outermost ring of pixels and
622
 *          (0.5 * fract * 255) for the second ring.  When blended over
623
 *          a second image, this
624
 *          (a) shrinks the visible image to make a clean overlap edge
625
 *              with an image below, and
626
 *          (b) softens the edges by weakening the aliasing there.
627
 *          Use l_setAlphaMaskBorder() to change these values.
628
 * </pre>
629
 */
630
PIX *
631
pixProjectivePtaWithAlpha(PIX       *pixs,
632
                          PTA       *ptad,
633
                          PTA       *ptas,
634
                          PIX       *pixg,
635
                          l_float32  fract,
636
                          l_int32    border)
637
0
{
638
0
l_int32  ws, hs, d;
639
0
PIX     *pixd, *pixb1, *pixb2, *pixg2, *pixga;
640
0
PTA     *ptad2, *ptas2;
641
642
0
    if (!pixs)
643
0
        return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
644
0
    pixGetDimensions(pixs, &ws, &hs, &d);
645
0
    if (d != 32 && pixGetColormap(pixs) == NULL)
646
0
        return (PIX *)ERROR_PTR("pixs not cmapped or 32 bpp", __func__, NULL);
647
0
    if (pixg && pixGetDepth(pixg) != 8) {
648
0
        L_WARNING("pixg not 8 bpp; using 'fract' transparent alpha\n",
649
0
                  __func__);
650
0
        pixg = NULL;
651
0
    }
652
0
    if (!pixg && (fract < 0.0 || fract > 1.0)) {
653
0
        L_WARNING("invalid fract; using 1.0 (fully transparent)\n", __func__);
654
0
        fract = 1.0;
655
0
    }
656
0
    if (!pixg && fract == 0.0)
657
0
        L_WARNING("fully opaque alpha; image will not be blended\n", __func__);
658
0
    if (!ptad)
659
0
        return (PIX *)ERROR_PTR("ptad not defined", __func__, NULL);
660
0
    if (!ptas)
661
0
        return (PIX *)ERROR_PTR("ptas not defined", __func__, NULL);
662
663
        /* Add border; the color doesn't matter */
664
0
    pixb1 = pixAddBorder(pixs, border, 0);
665
666
        /* Transform the ptr arrays to work on the bordered image */
667
0
    ptad2 = ptaTransform(ptad, border, border, 1.0, 1.0);
668
0
    ptas2 = ptaTransform(ptas, border, border, 1.0, 1.0);
669
670
        /* Do separate projective transform of rgb channels of pixs
671
         * and of pixg */
672
0
    pixd = pixProjectivePtaColor(pixb1, ptad2, ptas2, 0);
673
0
    if (!pixg) {
674
0
        pixg2 = pixCreate(ws, hs, 8);
675
0
        if (fract == 1.0)
676
0
            pixSetAll(pixg2);
677
0
        else
678
0
            pixSetAllArbitrary(pixg2, (l_int32)(255.0 * fract));
679
0
    } else {
680
0
        pixg2 = pixResizeToMatch(pixg, NULL, ws, hs);
681
0
    }
682
0
    if (ws > 10 && hs > 10) {  /* see note 7 */
683
0
        pixSetBorderRingVal(pixg2, 1,
684
0
                            (l_int32)(255.0 * fract * AlphaMaskBorderVals[0]));
685
0
        pixSetBorderRingVal(pixg2, 2,
686
0
                            (l_int32)(255.0 * fract * AlphaMaskBorderVals[1]));
687
688
0
    }
689
0
    pixb2 = pixAddBorder(pixg2, border, 0);  /* must be black border */
690
0
    pixga = pixProjectivePtaGray(pixb2, ptad2, ptas2, 0);
691
0
    pixSetRGBComponent(pixd, pixga, L_ALPHA_CHANNEL);
692
0
    pixSetSpp(pixd, 4);
693
694
0
    pixDestroy(&pixg2);
695
0
    pixDestroy(&pixb1);
696
0
    pixDestroy(&pixb2);
697
0
    pixDestroy(&pixga);
698
0
    ptaDestroy(&ptad2);
699
0
    ptaDestroy(&ptas2);
700
0
    return pixd;
701
0
}
702
703
704
/*-------------------------------------------------------------*
705
 *                Projective coordinate transformation         *
706
 *-------------------------------------------------------------*/
707
/*!
708
 * \brief   getProjectiveXformCoeffs()
709
 *
710
 * \param[in]    ptas   source 4 points; unprimed
711
 * \param[in]    ptad   transformed 4 points; primed
712
 * \param[out]   pvc    vector of coefficients of transform
713
 * \return  0 if OK; 1 on error
714
 *
715
 *  We have a set of 8 equations, describing the projective
716
 *  transformation that takes 4 points ptas into 4 other
717
 *  points ptad.  These equations are:
718
 *
719
 *          x1' = c[0]*x1 + c[1]*y1 + c[2]) / (c[6]*x1 + c[7]*y1 + 1
720
 *          y1' = c[3]*x1 + c[4]*y1 + c[5]) / (c[6]*x1 + c[7]*y1 + 1
721
 *          x2' = c[0]*x2 + c[1]*y2 + c[2]) / (c[6]*x2 + c[7]*y2 + 1
722
 *          y2' = c[3]*x2 + c[4]*y2 + c[5]) / (c[6]*x2 + c[7]*y2 + 1
723
 *          x3' = c[0]*x3 + c[1]*y3 + c[2]) / (c[6]*x3 + c[7]*y3 + 1
724
 *          y3' = c[3]*x3 + c[4]*y3 + c[5]) / (c[6]*x3 + c[7]*y3 + 1
725
 *          x4' = c[0]*x4 + c[1]*y4 + c[2]) / (c[6]*x4 + c[7]*y4 + 1
726
 *          y4' = c[3]*x4 + c[4]*y4 + c[5]) / (c[6]*x4 + c[7]*y4 + 1
727
 *
728
 *  Multiplying both sides of each eqn by the denominator, we get
729
 *
730
 *           AC = B
731
 *
732
 *  where B and C are column vectors
733
 *
734
 *         B = [ x1' y1' x2' y2' x3' y3' x4' y4' ]
735
 *         C = [ c[0] c[1] c[2] c[3] c[4] c[5] c[6] c[7] ]
736
 *
737
 *  and A is the 8x8 matrix
738
 *
739
 *             x1   y1     1     0   0    0   -x1*x1'  -y1*x1'
740
 *              0    0     0    x1   y1   1   -x1*y1'  -y1*y1'
741
 *             x2   y2     1     0   0    0   -x2*x2'  -y2*x2'
742
 *              0    0     0    x2   y2   1   -x2*y2'  -y2*y2'
743
 *             x3   y3     1     0   0    0   -x3*x3'  -y3*x3'
744
 *              0    0     0    x3   y3   1   -x3*y3'  -y3*y3'
745
 *             x4   y4     1     0   0    0   -x4*x4'  -y4*x4'
746
 *              0    0     0    x4   y4   1   -x4*y4'  -y4*y4'
747
 *
748
 *  These eight equations are solved here for the coefficients C.
749
 *
750
 *  These eight coefficients can then be used to find the mapping
751
 *  x,y) --> (x',y':
752
 *
753
 *           x' = c[0]x + c[1]y + c[2]) / (c[6]x + c[7]y + 1
754
 *           y' = c[3]x + c[4]y + c[5]) / (c[6]x + c[7]y + 1
755
 *
756
 *  that is implemented in projectiveXformSampled and
757
 *  projectiveXFormInterpolated.
758
 */
759
l_ok
760
getProjectiveXformCoeffs(PTA         *ptas,
761
                         PTA         *ptad,
762
                         l_float32  **pvc)
763
0
{
764
0
l_int32     i;
765
0
l_float32   x1, y1, x2, y2, x3, y3, x4, y4;
766
0
l_float32  *b;   /* rhs vector of primed coords X'; coeffs returned in *pvc */
767
0
l_float32  *a[8];  /* 8x8 matrix A  */
768
769
0
    if (!ptas)
770
0
        return ERROR_INT("ptas not defined", __func__, 1);
771
0
    if (!ptad)
772
0
        return ERROR_INT("ptad not defined", __func__, 1);
773
0
    if (!pvc)
774
0
        return ERROR_INT("&vc not defined", __func__, 1);
775
776
0
    b = (l_float32 *)LEPT_CALLOC(8, sizeof(l_float32));
777
0
    *pvc = b;
778
0
    ptaGetPt(ptas, 0, &x1, &y1);
779
0
    ptaGetPt(ptas, 1, &x2, &y2);
780
0
    ptaGetPt(ptas, 2, &x3, &y3);
781
0
    ptaGetPt(ptas, 3, &x4, &y4);
782
0
    ptaGetPt(ptad, 0, &b[0], &b[1]);
783
0
    ptaGetPt(ptad, 1, &b[2], &b[3]);
784
0
    ptaGetPt(ptad, 2, &b[4], &b[5]);
785
0
    ptaGetPt(ptad, 3, &b[6], &b[7]);
786
787
0
    for (i = 0; i < 8; i++)
788
0
        a[i] = (l_float32 *)LEPT_CALLOC(8, sizeof(l_float32));
789
0
    a[0][0] = x1;
790
0
    a[0][1] = y1;
791
0
    a[0][2] = 1.;
792
0
    a[0][6] = -x1 * b[0];
793
0
    a[0][7] = -y1 * b[0];
794
0
    a[1][3] = x1;
795
0
    a[1][4] = y1;
796
0
    a[1][5] = 1;
797
0
    a[1][6] = -x1 * b[1];
798
0
    a[1][7] = -y1 * b[1];
799
0
    a[2][0] = x2;
800
0
    a[2][1] = y2;
801
0
    a[2][2] = 1.;
802
0
    a[2][6] = -x2 * b[2];
803
0
    a[2][7] = -y2 * b[2];
804
0
    a[3][3] = x2;
805
0
    a[3][4] = y2;
806
0
    a[3][5] = 1;
807
0
    a[3][6] = -x2 * b[3];
808
0
    a[3][7] = -y2 * b[3];
809
0
    a[4][0] = x3;
810
0
    a[4][1] = y3;
811
0
    a[4][2] = 1.;
812
0
    a[4][6] = -x3 * b[4];
813
0
    a[4][7] = -y3 * b[4];
814
0
    a[5][3] = x3;
815
0
    a[5][4] = y3;
816
0
    a[5][5] = 1;
817
0
    a[5][6] = -x3 * b[5];
818
0
    a[5][7] = -y3 * b[5];
819
0
    a[6][0] = x4;
820
0
    a[6][1] = y4;
821
0
    a[6][2] = 1.;
822
0
    a[6][6] = -x4 * b[6];
823
0
    a[6][7] = -y4 * b[6];
824
0
    a[7][3] = x4;
825
0
    a[7][4] = y4;
826
0
    a[7][5] = 1;
827
0
    a[7][6] = -x4 * b[7];
828
0
    a[7][7] = -y4 * b[7];
829
830
0
    gaussjordan(a, b, 8);
831
832
0
    for (i = 0; i < 8; i++)
833
0
        LEPT_FREE(a[i]);
834
835
0
    return 0;
836
0
}
837
838
839
/*!
840
 * \brief   projectiveXformSampledPt()
841
 *
842
 * \param[in]    vc         vector of 8 coefficients
843
 * \param[in]    x, y       initial point
844
 * \param[out]   pxp, pyp   transformed point
845
 * \return  0 if OK; 1 on error
846
 *
847
 * <pre>
848
 * Notes:
849
 *      (1) This finds the nearest pixel coordinates of the transformed point.
850
 *      (2) It does not check ptrs for returned data!
851
 * </pre>
852
 */
853
l_ok
854
projectiveXformSampledPt(l_float32  *vc,
855
                         l_int32     x,
856
                         l_int32     y,
857
                         l_int32    *pxp,
858
                         l_int32    *pyp)
859
0
{
860
0
l_float32  factor;
861
0
l_float64  denom;
862
863
0
    if (!vc)
864
0
        return ERROR_INT("vc not defined", __func__, 1);
865
866
0
    if ((denom = vc[6] * x + vc[7] * y + 1.0f) == 0.0f)
867
0
        return ERROR_INT("denom = 0.0", __func__, 1);
868
0
    factor = 1.0f / denom;
869
0
    *pxp = (l_int32)(factor * (vc[0] * x + vc[1] * y + vc[2]) + 0.5f);
870
0
    *pyp = (l_int32)(factor * (vc[3] * x + vc[4] * y + vc[5]) + 0.5f);
871
0
    return 0;
872
0
}
873
874
875
/*!
876
 * \brief   projectiveXformPt()
877
 *
878
 * \param[in]    vc         vector of 8 coefficients
879
 * \param[in]    x, y       initial point
880
 * \param[out]   pxp, pyp   transformed point
881
 * \return  0 if OK; 1 on error
882
 *
883
 * <pre>
884
 * Notes:
885
 *      (1) This computes the floating point location of the transformed point.
886
 *      (2) It does not check ptrs for returned data!
887
 * </pre>
888
 */
889
l_ok
890
projectiveXformPt(l_float32  *vc,
891
                  l_int32     x,
892
                  l_int32     y,
893
                  l_float32  *pxp,
894
                  l_float32  *pyp)
895
0
{
896
0
l_float32  factor;
897
0
l_float64  denom;
898
899
0
    if (!vc)
900
0
        return ERROR_INT("vc not defined", __func__, 1);
901
902
0
    if ((denom = vc[6] * x + vc[7] * y + 1.0f) == 0.0f)
903
0
        return ERROR_INT("denom = 0.0", __func__, 1);
904
0
    factor = 1.0f / denom;
905
0
    *pxp = factor * (vc[0] * x + vc[1] * y + vc[2]);
906
0
    *pyp = factor * (vc[3] * x + vc[4] * y + vc[5]);
907
0
    return 0;
908
0
}