1/*M///////////////////////////////////////////////////////////////////////////////////////
2//
3//  IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
4//
5//  By downloading, copying, installing or using the software you agree to this license.
6//  If you do not agree to this license, do not download, install,
7//  copy or use the software.
8//
9//
10//                        Intel License Agreement
11//                For Open Source Computer Vision Library
12//
13// Copyright (C) 2000, Intel Corporation, all rights reserved.
14// Third party copyrights are property of their respective owners.
15//
16// Redistribution and use in source and binary forms, with or without modification,
17// are permitted provided that the following conditions are met:
18//
19//   * Redistribution's of source code must retain the above copyright notice,
20//     this list of conditions and the following disclaimer.
21//
22//   * Redistribution's in binary form must reproduce the above copyright notice,
23//     this list of conditions and the following disclaimer in the documentation
24//     and/or other materials provided with the distribution.
25//
26//   * The name of Intel Corporation may not be used to endorse or promote products
27//     derived from this software without specific prior written permission.
28//
29// This software is provided by the copyright holders and contributors "as is" and
30// any express or implied warranties, including, but not limited to, the implied
31// warranties of merchantability and fitness for a particular purpose are disclaimed.
32// In no event shall the Intel Corporation or contributors be liable for any direct,
33// indirect, incidental, special, exemplary, or consequential damages
34// (including, but not limited to, procurement of substitute goods or services;
35// loss of use, data, or profits; or business interruption) however caused
36// and on any theory of liability, whether in contract, strict liability,
37// or tort (including negligence or otherwise) arising in any way out of
38// the use of this software, even if advised of the possibility of such damage.
39//
40//M*/
41#include "_cv.h"
42
43#define ICV_DIST_SHIFT  16
44#define ICV_INIT_DIST0  (INT_MAX >> 2)
45
46static CvStatus
47icvInitTopBottom( int* temp, int tempstep, CvSize size, int border )
48{
49    int i, j;
50    for( i = 0; i < border; i++ )
51    {
52        int* ttop = (int*)(temp + i*tempstep);
53        int* tbottom = (int*)(temp + (size.height + border*2 - i - 1)*tempstep);
54
55        for( j = 0; j < size.width + border*2; j++ )
56        {
57            ttop[j] = ICV_INIT_DIST0;
58            tbottom[j] = ICV_INIT_DIST0;
59        }
60    }
61
62    return CV_OK;
63}
64
65
66static CvStatus CV_STDCALL
67icvDistanceTransform_3x3_C1R( const uchar* src, int srcstep, int* temp,
68        int step, float* dist, int dststep, CvSize size, const float* metrics )
69{
70    const int BORDER = 1;
71    int i, j;
72    const int HV_DIST = CV_FLT_TO_FIX( metrics[0], ICV_DIST_SHIFT );
73    const int DIAG_DIST = CV_FLT_TO_FIX( metrics[1], ICV_DIST_SHIFT );
74    const float scale = 1.f/(1 << ICV_DIST_SHIFT);
75
76    srcstep /= sizeof(src[0]);
77    step /= sizeof(temp[0]);
78    dststep /= sizeof(dist[0]);
79
80    icvInitTopBottom( temp, step, size, BORDER );
81
82    // forward pass
83    for( i = 0; i < size.height; i++ )
84    {
85        const uchar* s = src + i*srcstep;
86        int* tmp = (int*)(temp + (i+BORDER)*step) + BORDER;
87
88        for( j = 0; j < BORDER; j++ )
89            tmp[-j-1] = tmp[size.width + j] = ICV_INIT_DIST0;
90
91        for( j = 0; j < size.width; j++ )
92        {
93            if( !s[j] )
94                tmp[j] = 0;
95            else
96            {
97                int t0 = tmp[j-step-1] + DIAG_DIST;
98                int t = tmp[j-step] + HV_DIST;
99                if( t0 > t ) t0 = t;
100                t = tmp[j-step+1] + DIAG_DIST;
101                if( t0 > t ) t0 = t;
102                t = tmp[j-1] + HV_DIST;
103                if( t0 > t ) t0 = t;
104                tmp[j] = t0;
105            }
106        }
107    }
108
109    // backward pass
110    for( i = size.height - 1; i >= 0; i-- )
111    {
112        float* d = (float*)(dist + i*dststep);
113        int* tmp = (int*)(temp + (i+BORDER)*step) + BORDER;
114
115        for( j = size.width - 1; j >= 0; j-- )
116        {
117            int t0 = tmp[j];
118            if( t0 > HV_DIST )
119            {
120                int t = tmp[j+step+1] + DIAG_DIST;
121                if( t0 > t ) t0 = t;
122                t = tmp[j+step] + HV_DIST;
123                if( t0 > t ) t0 = t;
124                t = tmp[j+step-1] + DIAG_DIST;
125                if( t0 > t ) t0 = t;
126                t = tmp[j+1] + HV_DIST;
127                if( t0 > t ) t0 = t;
128                tmp[j] = t0;
129            }
130            d[j] = (float)(t0 * scale);
131        }
132    }
133
134    return CV_OK;
135}
136
137
138static CvStatus CV_STDCALL
139icvDistanceTransform_5x5_C1R( const uchar* src, int srcstep, int* temp,
140        int step, float* dist, int dststep, CvSize size, const float* metrics )
141{
142    const int BORDER = 2;
143    int i, j;
144    const int HV_DIST = CV_FLT_TO_FIX( metrics[0], ICV_DIST_SHIFT );
145    const int DIAG_DIST = CV_FLT_TO_FIX( metrics[1], ICV_DIST_SHIFT );
146    const int LONG_DIST = CV_FLT_TO_FIX( metrics[2], ICV_DIST_SHIFT );
147    const float scale = 1.f/(1 << ICV_DIST_SHIFT);
148
149    srcstep /= sizeof(src[0]);
150    step /= sizeof(temp[0]);
151    dststep /= sizeof(dist[0]);
152
153    icvInitTopBottom( temp, step, size, BORDER );
154
155    // forward pass
156    for( i = 0; i < size.height; i++ )
157    {
158        const uchar* s = src + i*srcstep;
159        int* tmp = (int*)(temp + (i+BORDER)*step) + BORDER;
160
161        for( j = 0; j < BORDER; j++ )
162            tmp[-j-1] = tmp[size.width + j] = ICV_INIT_DIST0;
163
164        for( j = 0; j < size.width; j++ )
165        {
166            if( !s[j] )
167                tmp[j] = 0;
168            else
169            {
170                int t0 = tmp[j-step*2-1] + LONG_DIST;
171                int t = tmp[j-step*2+1] + LONG_DIST;
172                if( t0 > t ) t0 = t;
173                t = tmp[j-step-2] + LONG_DIST;
174                if( t0 > t ) t0 = t;
175                t = tmp[j-step-1] + DIAG_DIST;
176                if( t0 > t ) t0 = t;
177                t = tmp[j-step] + HV_DIST;
178                if( t0 > t ) t0 = t;
179                t = tmp[j-step+1] + DIAG_DIST;
180                if( t0 > t ) t0 = t;
181                t = tmp[j-step+2] + LONG_DIST;
182                if( t0 > t ) t0 = t;
183                t = tmp[j-1] + HV_DIST;
184                if( t0 > t ) t0 = t;
185                tmp[j] = t0;
186            }
187        }
188    }
189
190    // backward pass
191    for( i = size.height - 1; i >= 0; i-- )
192    {
193        float* d = (float*)(dist + i*dststep);
194        int* tmp = (int*)(temp + (i+BORDER)*step) + BORDER;
195
196        for( j = size.width - 1; j >= 0; j-- )
197        {
198            int t0 = tmp[j];
199            if( t0 > HV_DIST )
200            {
201                int t = tmp[j+step*2+1] + LONG_DIST;
202                if( t0 > t ) t0 = t;
203                t = tmp[j+step*2-1] + LONG_DIST;
204                if( t0 > t ) t0 = t;
205                t = tmp[j+step+2] + LONG_DIST;
206                if( t0 > t ) t0 = t;
207                t = tmp[j+step+1] + DIAG_DIST;
208                if( t0 > t ) t0 = t;
209                t = tmp[j+step] + HV_DIST;
210                if( t0 > t ) t0 = t;
211                t = tmp[j+step-1] + DIAG_DIST;
212                if( t0 > t ) t0 = t;
213                t = tmp[j+step-2] + LONG_DIST;
214                if( t0 > t ) t0 = t;
215                t = tmp[j+1] + HV_DIST;
216                if( t0 > t ) t0 = t;
217                tmp[j] = t0;
218            }
219            d[j] = (float)(t0 * scale);
220        }
221    }
222
223    return CV_OK;
224}
225
226
227static CvStatus CV_STDCALL
228icvDistanceTransformEx_5x5_C1R( const uchar* src, int srcstep, int* temp,
229                int step, float* dist, int dststep, int* labels, int lstep,
230                CvSize size, const float* metrics )
231{
232    const int BORDER = 2;
233
234    int i, j;
235    const int HV_DIST = CV_FLT_TO_FIX( metrics[0], ICV_DIST_SHIFT );
236    const int DIAG_DIST = CV_FLT_TO_FIX( metrics[1], ICV_DIST_SHIFT );
237    const int LONG_DIST = CV_FLT_TO_FIX( metrics[2], ICV_DIST_SHIFT );
238    const float scale = 1.f/(1 << ICV_DIST_SHIFT);
239
240    srcstep /= sizeof(src[0]);
241    step /= sizeof(temp[0]);
242    dststep /= sizeof(dist[0]);
243    lstep /= sizeof(labels[0]);
244
245    icvInitTopBottom( temp, step, size, BORDER );
246
247    // forward pass
248    for( i = 0; i < size.height; i++ )
249    {
250        const uchar* s = src + i*srcstep;
251        int* tmp = (int*)(temp + (i+BORDER)*step) + BORDER;
252        int* lls = (int*)(labels + i*lstep);
253
254        for( j = 0; j < BORDER; j++ )
255            tmp[-j-1] = tmp[size.width + j] = ICV_INIT_DIST0;
256
257        for( j = 0; j < size.width; j++ )
258        {
259            if( !s[j] )
260            {
261                tmp[j] = 0;
262                //assert( lls[j] != 0 );
263            }
264            else
265            {
266                int t0 = ICV_INIT_DIST0, t;
267                int l0 = 0;
268
269                t = tmp[j-step*2-1] + LONG_DIST;
270                if( t0 > t )
271                {
272                    t0 = t;
273                    l0 = lls[j-lstep*2-1];
274                }
275                t = tmp[j-step*2+1] + LONG_DIST;
276                if( t0 > t )
277                {
278                    t0 = t;
279                    l0 = lls[j-lstep*2+1];
280                }
281                t = tmp[j-step-2] + LONG_DIST;
282                if( t0 > t )
283                {
284                    t0 = t;
285                    l0 = lls[j-lstep-2];
286                }
287                t = tmp[j-step-1] + DIAG_DIST;
288                if( t0 > t )
289                {
290                    t0 = t;
291                    l0 = lls[j-lstep-1];
292                }
293                t = tmp[j-step] + HV_DIST;
294                if( t0 > t )
295                {
296                    t0 = t;
297                    l0 = lls[j-lstep];
298                }
299                t = tmp[j-step+1] + DIAG_DIST;
300                if( t0 > t )
301                {
302                    t0 = t;
303                    l0 = lls[j-lstep+1];
304                }
305                t = tmp[j-step+2] + LONG_DIST;
306                if( t0 > t )
307                {
308                    t0 = t;
309                    l0 = lls[j-lstep+2];
310                }
311                t = tmp[j-1] + HV_DIST;
312                if( t0 > t )
313                {
314                    t0 = t;
315                    l0 = lls[j-1];
316                }
317
318                tmp[j] = t0;
319                lls[j] = l0;
320            }
321        }
322    }
323
324    // backward pass
325    for( i = size.height - 1; i >= 0; i-- )
326    {
327        float* d = (float*)(dist + i*dststep);
328        int* tmp = (int*)(temp + (i+BORDER)*step) + BORDER;
329        int* lls = (int*)(labels + i*lstep);
330
331        for( j = size.width - 1; j >= 0; j-- )
332        {
333            int t0 = tmp[j];
334            int l0 = lls[j];
335            if( t0 > HV_DIST )
336            {
337                int t = tmp[j+step*2+1] + LONG_DIST;
338                if( t0 > t )
339                {
340                    t0 = t;
341                    l0 = lls[j+lstep*2+1];
342                }
343                t = tmp[j+step*2-1] + LONG_DIST;
344                if( t0 > t )
345                {
346                    t0 = t;
347                    l0 = lls[j+lstep*2-1];
348                }
349                t = tmp[j+step+2] + LONG_DIST;
350                if( t0 > t )
351                {
352                    t0 = t;
353                    l0 = lls[j+lstep+2];
354                }
355                t = tmp[j+step+1] + DIAG_DIST;
356                if( t0 > t )
357                {
358                    t0 = t;
359                    l0 = lls[j+lstep+1];
360                }
361                t = tmp[j+step] + HV_DIST;
362                if( t0 > t )
363                {
364                    t0 = t;
365                    l0 = lls[j+lstep];
366                }
367                t = tmp[j+step-1] + DIAG_DIST;
368                if( t0 > t )
369                {
370                    t0 = t;
371                    l0 = lls[j+lstep-1];
372                }
373                t = tmp[j+step-2] + LONG_DIST;
374                if( t0 > t )
375                {
376                    t0 = t;
377                    l0 = lls[j+lstep-2];
378                }
379                t = tmp[j+1] + HV_DIST;
380                if( t0 > t )
381                {
382                    t0 = t;
383                    l0 = lls[j+1];
384                }
385                tmp[j] = t0;
386                lls[j] = l0;
387            }
388            d[j] = (float)(t0 * scale);
389        }
390    }
391
392    return CV_OK;
393}
394
395
396static CvStatus
397icvGetDistanceTransformMask( int maskType, float *metrics )
398{
399    if( !metrics )
400        return CV_NULLPTR_ERR;
401
402    switch (maskType)
403    {
404    case 30:
405        metrics[0] = 1.0f;
406        metrics[1] = 1.0f;
407        break;
408
409    case 31:
410        metrics[0] = 1.0f;
411        metrics[1] = 2.0f;
412        break;
413
414    case 32:
415        metrics[0] = 0.955f;
416        metrics[1] = 1.3693f;
417        break;
418
419    case 50:
420        metrics[0] = 1.0f;
421        metrics[1] = 1.0f;
422        metrics[2] = 2.0f;
423        break;
424
425    case 51:
426        metrics[0] = 1.0f;
427        metrics[1] = 2.0f;
428        metrics[2] = 3.0f;
429        break;
430
431    case 52:
432        metrics[0] = 1.0f;
433        metrics[1] = 1.4f;
434        metrics[2] = 2.1969f;
435        break;
436    default:
437        return CV_BADRANGE_ERR;
438    }
439
440    return CV_OK;
441}
442
443
444static void
445icvTrueDistTrans( const CvMat* src, CvMat* dst )
446{
447    CvMat* buffer = 0;
448
449    CV_FUNCNAME( "cvDistTransform2" );
450
451    __BEGIN__;
452
453    int i, m, n;
454    int sstep, dstep;
455    const float inf = 1e6f;
456    int thread_count = cvGetNumThreads();
457    int pass1_sz, pass2_sz;
458
459    if( !CV_ARE_SIZES_EQ( src, dst ))
460        CV_ERROR( CV_StsUnmatchedSizes, "" );
461
462    if( CV_MAT_TYPE(src->type) != CV_8UC1 ||
463        CV_MAT_TYPE(dst->type) != CV_32FC1 )
464        CV_ERROR( CV_StsUnsupportedFormat,
465        "The input image must have 8uC1 type and the output one must have 32fC1 type" );
466
467    m = src->rows;
468    n = src->cols;
469
470    // (see stage 1 below):
471    // sqr_tab: 2*m, sat_tab: 3*m + 1, d: m*thread_count,
472    pass1_sz = src->rows*(5 + thread_count) + 1;
473    // (see stage 2):
474    // sqr_tab & inv_tab: n each; f & v: n*thread_count each; z: (n+1)*thread_count
475    pass2_sz = src->cols*(2 + thread_count*3) + thread_count;
476    CV_CALL( buffer = cvCreateMat( 1, MAX(pass1_sz, pass2_sz), CV_32FC1 ));
477
478    sstep = src->step;
479    dstep = dst->step / sizeof(float);
480
481    // stage 1: compute 1d distance transform of each column
482    {
483    float* sqr_tab = buffer->data.fl;
484    int* sat_tab = (int*)(sqr_tab + m*2);
485    const int shift = m*2;
486
487    for( i = 0; i < m; i++ )
488        sqr_tab[i] = (float)(i*i);
489    for( i = m; i < m*2; i++ )
490        sqr_tab[i] = inf;
491    for( i = 0; i < shift; i++ )
492        sat_tab[i] = 0;
493    for( ; i <= m*3; i++ )
494        sat_tab[i] = i - shift;
495
496#ifdef _OPENMP
497    #pragma omp parallel for num_threads(thread_count)
498#endif
499    for( i = 0; i < n; i++ )
500    {
501        const uchar* sptr = src->data.ptr + i + (m-1)*sstep;
502        float* dptr = dst->data.fl + i;
503        int* d = (int*)(sat_tab + m*3+1+m*cvGetThreadNum());
504        int j, dist = m-1;
505
506        for( j = m-1; j >= 0; j--, sptr -= sstep )
507        {
508            dist = (dist + 1) & (sptr[0] == 0 ? 0 : -1);
509            d[j] = dist;
510        }
511
512        dist = m-1;
513        for( j = 0; j < m; j++, dptr += dstep )
514        {
515            dist = dist + 1 - sat_tab[dist + 1 - d[j] + shift];
516            d[j] = dist;
517            dptr[0] = sqr_tab[dist];
518        }
519    }
520    }
521
522    // stage 2: compute modified distance transform for each row
523    {
524    float* inv_tab = buffer->data.fl;
525    float* sqr_tab = inv_tab + n;
526
527    inv_tab[0] = sqr_tab[0] = 0.f;
528    for( i = 1; i < n; i++ )
529    {
530        inv_tab[i] = (float)(0.5/i);
531        sqr_tab[i] = (float)(i*i);
532    }
533
534#ifdef _OPENMP
535    #pragma omp parallel for num_threads(thread_count) schedule(dynamic)
536#endif
537    for( i = 0; i < m; i++ )
538    {
539        float* d = (float*)(dst->data.ptr + i*dst->step);
540        float* f = sqr_tab + n + (n*3+1)*cvGetThreadNum();
541        float* z = f + n;
542        int* v = (int*)(z + n + 1);
543        int p, q, k;
544
545        v[0] = 0;
546        z[0] = -inf;
547        z[1] = inf;
548        f[0] = d[0];
549
550        for( q = 1, k = 0; q < n; q++ )
551        {
552            float fq = d[q];
553            f[q] = fq;
554
555            for(;;k--)
556            {
557                p = v[k];
558                float s = (fq + sqr_tab[q] - d[p] - sqr_tab[p])*inv_tab[q - p];
559                if( s > z[k] )
560                {
561                    k++;
562                    v[k] = q;
563                    z[k] = s;
564                    z[k+1] = inf;
565                    break;
566                }
567            }
568        }
569
570        for( q = 0, k = 0; q < n; q++ )
571        {
572            while( z[k+1] < q )
573                k++;
574            p = v[k];
575            d[q] = sqr_tab[abs(q - p)] + f[p];
576        }
577    }
578    }
579
580    cvPow( dst, dst, 0.5 );
581
582    __END__;
583
584    cvReleaseMat( &buffer );
585}
586
587
588/*********************************** IPP functions *********************************/
589
590icvDistanceTransform_3x3_8u32f_C1R_t icvDistanceTransform_3x3_8u32f_C1R_p = 0;
591icvDistanceTransform_5x5_8u32f_C1R_t icvDistanceTransform_5x5_8u32f_C1R_p = 0;
592icvDistanceTransform_3x3_8u_C1IR_t icvDistanceTransform_3x3_8u_C1IR_p = 0;
593icvDistanceTransform_3x3_8u_C1R_t icvDistanceTransform_3x3_8u_C1R_p = 0;
594
595typedef CvStatus (CV_STDCALL * CvIPPDistTransFunc)( const uchar* src, int srcstep,
596                                                    void* dst, int dststep,
597                                                    CvSize size, const void* metrics );
598
599typedef CvStatus (CV_STDCALL * CvIPPDistTransFunc2)( uchar* src, int srcstep,
600                                                     CvSize size, const int* metrics );
601
602/***********************************************************************************/
603
604typedef CvStatus (CV_STDCALL * CvDistTransFunc)( const uchar* src, int srcstep,
605                                                 int* temp, int tempstep,
606                                                 float* dst, int dststep,
607                                                 CvSize size, const float* metrics );
608
609
610/****************************************************************************************\
611 User-contributed code:
612
613 Non-inplace and Inplace 8u->8u Distance Transform for CityBlock (a.k.a. L1) metric
614 (C) 2006 by Jay Stavinzky.
615\****************************************************************************************/
616
617//BEGIN ATS ADDITION
618/* 8-bit grayscale distance transform function */
619static void
620icvDistanceATS_L1_8u( const CvMat* src, CvMat* dst )
621{
622    CV_FUNCNAME( "cvDistanceATS" );
623
624    __BEGIN__;
625
626    int width = src->cols, height = src->rows;
627
628    int a;
629    uchar lut[256];
630    int x, y;
631
632    const uchar *sbase = src->data.ptr;
633    uchar *dbase = dst->data.ptr;
634    int srcstep = src->step;
635    int dststep = dst->step;
636
637    CV_ASSERT( CV_IS_MASK_ARR( src ) && CV_MAT_TYPE( dst->type ) == CV_8UC1 );
638    CV_ASSERT( CV_ARE_SIZES_EQ( src, dst ));
639
640    ////////////////////// forward scan ////////////////////////
641    for( x = 0; x < 256; x++ )
642        lut[x] = CV_CAST_8U(x+1);
643
644    //init first pixel to max (we're going to be skipping it)
645    dbase[0] = (uchar)(sbase[0] == 0 ? 0 : 255);
646
647    //first row (scan west only, skip first pixel)
648    for( x = 1; x < width; x++ )
649        dbase[x] = (uchar)(sbase[x] == 0 ? 0 : lut[dbase[x-1]]);
650
651    for( y = 1; y < height; y++ )
652    {
653        sbase += srcstep;
654        dbase += dststep;
655
656        //for left edge, scan north only
657        a = sbase[0] == 0 ? 0 : lut[dbase[-dststep]];
658        dbase[0] = (uchar)a;
659
660        for( x = 1; x < width; x++ )
661        {
662            a = sbase[x] == 0 ? 0 : lut[MIN(a, dbase[x - dststep])];
663            dbase[x] = (uchar)a;
664        }
665    }
666
667    ////////////////////// backward scan ///////////////////////
668
669    a = dbase[width-1];
670
671    // do last row east pixel scan here (skip bottom right pixel)
672    for( x = width - 2; x >= 0; x-- )
673    {
674        a = lut[a];
675        dbase[x] = (uchar)(CV_CALC_MIN_8U(a, dbase[x]));
676    }
677
678    // right edge is the only error case
679    for( y = height - 2; y >= 0; y-- )
680    {
681        dbase -= dststep;
682
683        // do right edge
684        a = lut[dbase[width-1+dststep]];
685        dbase[width-1] = (uchar)(MIN(a, dbase[width-1]));
686
687        for( x = width - 2; x >= 0; x-- )
688        {
689            int b = dbase[x+dststep];
690            a = lut[MIN(a, b)];
691            dbase[x] = (uchar)(MIN(a, dbase[x]));
692        }
693    }
694
695    __END__;
696}
697//END ATS ADDITION
698
699
700/* Wrapper function for distance transform group */
701CV_IMPL void
702cvDistTransform( const void* srcarr, void* dstarr,
703                 int distType, int maskSize,
704                 const float *mask,
705                 void* labelsarr )
706{
707    CvMat* temp = 0;
708    CvMat* src_copy = 0;
709    CvMemStorage* st = 0;
710
711    CV_FUNCNAME( "cvDistTransform" );
712
713    __BEGIN__;
714
715    float _mask[5] = {0};
716    int _imask[3];
717    CvMat srcstub, *src = (CvMat*)srcarr;
718    CvMat dststub, *dst = (CvMat*)dstarr;
719    CvMat lstub, *labels = (CvMat*)labelsarr;
720    CvSize size;
721    CvIPPDistTransFunc ipp_func = 0;
722    CvIPPDistTransFunc2 ipp_inp_func = 0;
723
724    CV_CALL( src = cvGetMat( src, &srcstub ));
725    CV_CALL( dst = cvGetMat( dst, &dststub ));
726
727    if( !CV_IS_MASK_ARR( src ) || (CV_MAT_TYPE( dst->type ) != CV_32FC1 &&
728        (CV_MAT_TYPE(dst->type) != CV_8UC1 || distType != CV_DIST_L1 || labels)) )
729        CV_ERROR( CV_StsUnsupportedFormat,
730        "source image must be 8uC1 and the distance map must be 32fC1 "
731        "(or 8uC1 in case of simple L1 distance transform)" );
732
733    if( !CV_ARE_SIZES_EQ( src, dst ))
734        CV_ERROR( CV_StsUnmatchedSizes, "the source and the destination images must be of the same size" );
735
736    if( maskSize != CV_DIST_MASK_3 && maskSize != CV_DIST_MASK_5 && maskSize != CV_DIST_MASK_PRECISE )
737        CV_ERROR( CV_StsBadSize, "Mask size should be 3 or 5 or 0 (presize)" );
738
739    if( distType == CV_DIST_C || distType == CV_DIST_L1 )
740        maskSize = !labels ? CV_DIST_MASK_3 : CV_DIST_MASK_5;
741    else if( distType == CV_DIST_L2 && labels )
742        maskSize = CV_DIST_MASK_5;
743
744    if( maskSize == CV_DIST_MASK_PRECISE )
745    {
746        CV_CALL( icvTrueDistTrans( src, dst ));
747        EXIT;
748    }
749
750    if( labels )
751    {
752        CV_CALL( labels = cvGetMat( labels, &lstub ));
753        if( CV_MAT_TYPE( labels->type ) != CV_32SC1 )
754            CV_ERROR( CV_StsUnsupportedFormat, "the output array of labels must be 32sC1" );
755
756        if( !CV_ARE_SIZES_EQ( labels, dst ))
757            CV_ERROR( CV_StsUnmatchedSizes, "the array of labels has a different size" );
758
759        if( maskSize == CV_DIST_MASK_3 )
760            CV_ERROR( CV_StsNotImplemented,
761            "3x3 mask can not be used for \"labeled\" distance transform. Use 5x5 mask" );
762    }
763
764    if( distType == CV_DIST_C || distType == CV_DIST_L1 || distType == CV_DIST_L2 )
765    {
766        icvGetDistanceTransformMask( (distType == CV_DIST_C ? 0 :
767            distType == CV_DIST_L1 ? 1 : 2) + maskSize*10, _mask );
768    }
769    else if( distType == CV_DIST_USER )
770    {
771        if( !mask )
772            CV_ERROR( CV_StsNullPtr, "" );
773
774        memcpy( _mask, mask, (maskSize/2 + 1)*sizeof(float));
775    }
776
777    if( !labels )
778    {
779        if( CV_MAT_TYPE(dst->type) == CV_32FC1 )
780            ipp_func = (CvIPPDistTransFunc)(maskSize == CV_DIST_MASK_3 ?
781                icvDistanceTransform_3x3_8u32f_C1R_p : icvDistanceTransform_5x5_8u32f_C1R_p);
782        else if( src->data.ptr != dst->data.ptr )
783            ipp_func = (CvIPPDistTransFunc)icvDistanceTransform_3x3_8u_C1R_p;
784        else
785            ipp_inp_func = icvDistanceTransform_3x3_8u_C1IR_p;
786    }
787
788    size = cvGetMatSize(src);
789
790    if( (ipp_func || ipp_inp_func) && src->cols >= 4 && src->rows >= 2 )
791    {
792        _imask[0] = cvRound(_mask[0]);
793        _imask[1] = cvRound(_mask[1]);
794        _imask[2] = cvRound(_mask[2]);
795
796        if( ipp_func )
797        {
798            IPPI_CALL( ipp_func( src->data.ptr, src->step,
799                    dst->data.fl, dst->step, size,
800                    CV_MAT_TYPE(dst->type) == CV_8UC1 ?
801                    (void*)_imask : (void*)_mask ));
802        }
803        else
804        {
805            IPPI_CALL( ipp_inp_func( src->data.ptr, src->step, size, _imask ));
806        }
807    }
808    else if( CV_MAT_TYPE(dst->type) == CV_8UC1 )
809    {
810        CV_CALL( icvDistanceATS_L1_8u( src, dst ));
811    }
812    else
813    {
814        int border = maskSize == CV_DIST_MASK_3 ? 1 : 2;
815        CV_CALL( temp = cvCreateMat( size.height + border*2, size.width + border*2, CV_32SC1 ));
816
817        if( !labels )
818        {
819            CvDistTransFunc func = maskSize == CV_DIST_MASK_3 ?
820                icvDistanceTransform_3x3_C1R :
821                icvDistanceTransform_5x5_C1R;
822
823            func( src->data.ptr, src->step, temp->data.i, temp->step,
824                  dst->data.fl, dst->step, size, _mask );
825        }
826        else
827        {
828            CvSeq *contours = 0;
829            CvPoint top_left = {0,0}, bottom_right = {size.width-1,size.height-1};
830            int label;
831
832            CV_CALL( st = cvCreateMemStorage() );
833            CV_CALL( src_copy = cvCreateMat( size.height, size.width, src->type ));
834            cvCmpS( src, 0, src_copy, CV_CMP_EQ );
835            cvFindContours( src_copy, st, &contours, sizeof(CvContour),
836                            CV_RETR_CCOMP, CV_CHAIN_APPROX_SIMPLE );
837            cvZero( labels );
838            for( label = 1; contours != 0; contours = contours->h_next, label++ )
839            {
840                CvScalar area_color = cvScalarAll(label);
841                cvDrawContours( labels, contours, area_color, area_color, -255, -1, 8 );
842            }
843
844            cvCopy( src, src_copy );
845            cvRectangle( src_copy, top_left, bottom_right, cvScalarAll(255), 1, 8 );
846
847            icvDistanceTransformEx_5x5_C1R( src_copy->data.ptr, src_copy->step, temp->data.i, temp->step,
848                        dst->data.fl, dst->step, labels->data.i, labels->step, size, _mask );
849        }
850    }
851
852    __END__;
853
854    cvReleaseMat( &temp );
855    cvReleaseMat( &src_copy );
856    cvReleaseMemStorage( &st );
857}
858
859/* End of file. */
860