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
42/****************************************************************************************\
43
44      Calculation of a texture descriptors from GLCM (Grey Level Co-occurrence Matrix'es)
45      The code was submitted by Daniel Eaton [danieljameseaton@yahoo.com]
46
47\****************************************************************************************/
48
49#include "_cvaux.h"
50
51#include <math.h>
52#include <assert.h>
53
54#define CV_MAX_NUM_GREY_LEVELS_8U  256
55
56struct CvGLCM
57{
58    int matrixSideLength;
59    int numMatrices;
60    double*** matrices;
61
62    int  numLookupTableElements;
63    int  forwardLookupTable[CV_MAX_NUM_GREY_LEVELS_8U];
64    int  reverseLookupTable[CV_MAX_NUM_GREY_LEVELS_8U];
65
66    double** descriptors;
67    int numDescriptors;
68    int descriptorOptimizationType;
69    int optimizationType;
70};
71
72
73static void icvCreateGLCM_LookupTable_8u_C1R( const uchar* srcImageData, int srcImageStep,
74                                             CvSize srcImageSize, CvGLCM* destGLCM,
75                                             int* steps, int numSteps, int* memorySteps );
76
77static void
78icvCreateGLCMDescriptors_AllowDoubleNest( CvGLCM* destGLCM, int matrixIndex );
79
80
81CV_IMPL CvGLCM*
82cvCreateGLCM( const IplImage* srcImage,
83              int stepMagnitude,
84              const int* srcStepDirections,/* should be static array..
85                                          or if not the user should handle de-allocation */
86              int numStepDirections,
87              int optimizationType )
88{
89    static const int defaultStepDirections[] = { 0,1, -1,1, -1,0, -1,-1 };
90
91    int* memorySteps = 0;
92    CvGLCM* newGLCM = 0;
93    int* stepDirections = 0;
94
95    CV_FUNCNAME( "cvCreateGLCM" );
96
97    __BEGIN__;
98
99    uchar* srcImageData = 0;
100    CvSize srcImageSize;
101    int srcImageStep;
102    int stepLoop;
103    const int maxNumGreyLevels8u = CV_MAX_NUM_GREY_LEVELS_8U;
104
105    if( !srcImage )
106        CV_ERROR( CV_StsNullPtr, "" );
107
108    if( srcImage->nChannels != 1 )
109        CV_ERROR( CV_BadNumChannels, "Number of channels must be 1");
110
111    if( srcImage->depth != IPL_DEPTH_8U )
112        CV_ERROR( CV_BadDepth, "Depth must be equal IPL_DEPTH_8U");
113
114    // no Directions provided, use the default ones - 0 deg, 45, 90, 135
115    if( !srcStepDirections )
116    {
117        srcStepDirections = defaultStepDirections;
118    }
119
120    CV_CALL( stepDirections = (int*)cvAlloc( numStepDirections*2*sizeof(stepDirections[0])));
121    memcpy( stepDirections, srcStepDirections, numStepDirections*2*sizeof(stepDirections[0]));
122
123    cvGetImageRawData( srcImage, &srcImageData, &srcImageStep, &srcImageSize );
124
125    // roll together Directions and magnitudes together with knowledge of image (step)
126    CV_CALL( memorySteps = (int*)cvAlloc( numStepDirections*sizeof(memorySteps[0])));
127
128    for( stepLoop = 0; stepLoop < numStepDirections; stepLoop++ )
129    {
130        stepDirections[stepLoop*2 + 0] *= stepMagnitude;
131        stepDirections[stepLoop*2 + 1] *= stepMagnitude;
132
133        memorySteps[stepLoop] = stepDirections[stepLoop*2 + 0]*srcImageStep +
134                                stepDirections[stepLoop*2 + 1];
135    }
136
137    CV_CALL( newGLCM = (CvGLCM*)cvAlloc(sizeof(newGLCM)));
138    memset( newGLCM, 0, sizeof(newGLCM) );
139
140    newGLCM->matrices = 0;
141    newGLCM->numMatrices = numStepDirections;
142    newGLCM->optimizationType = optimizationType;
143
144    if( optimizationType <= CV_GLCM_OPTIMIZATION_LUT )
145    {
146        int lookupTableLoop, imageColLoop, imageRowLoop, lineOffset = 0;
147
148        // if optimization type is set to lut, then make one for the image
149        if( optimizationType == CV_GLCM_OPTIMIZATION_LUT )
150        {
151            for( imageRowLoop = 0; imageRowLoop < srcImageSize.height;
152                                   imageRowLoop++, lineOffset += srcImageStep )
153            {
154                for( imageColLoop = 0; imageColLoop < srcImageSize.width; imageColLoop++ )
155                {
156                    newGLCM->forwardLookupTable[srcImageData[lineOffset+imageColLoop]]=1;
157                }
158            }
159
160            newGLCM->numLookupTableElements = 0;
161
162            for( lookupTableLoop = 0; lookupTableLoop < maxNumGreyLevels8u; lookupTableLoop++ )
163            {
164                if( newGLCM->forwardLookupTable[ lookupTableLoop ] != 0 )
165                {
166                    newGLCM->forwardLookupTable[ lookupTableLoop ] =
167                        newGLCM->numLookupTableElements;
168                    newGLCM->reverseLookupTable[ newGLCM->numLookupTableElements ] =
169                        lookupTableLoop;
170
171                    newGLCM->numLookupTableElements++;
172                }
173            }
174        }
175        // otherwise make a "LUT" which contains all the gray-levels (for code-reuse)
176        else if( optimizationType == CV_GLCM_OPTIMIZATION_NONE )
177        {
178            for( lookupTableLoop = 0; lookupTableLoop <maxNumGreyLevels8u; lookupTableLoop++ )
179            {
180                newGLCM->forwardLookupTable[ lookupTableLoop ] = lookupTableLoop;
181                newGLCM->reverseLookupTable[ lookupTableLoop ] = lookupTableLoop;
182            }
183            newGLCM->numLookupTableElements = maxNumGreyLevels8u;
184        }
185
186        newGLCM->matrixSideLength = newGLCM->numLookupTableElements;
187        icvCreateGLCM_LookupTable_8u_C1R( srcImageData, srcImageStep, srcImageSize,
188                                          newGLCM, stepDirections,
189                                          numStepDirections, memorySteps );
190    }
191    else if( optimizationType == CV_GLCM_OPTIMIZATION_HISTOGRAM )
192    {
193        CV_ERROR( CV_StsBadFlag, "Histogram-based method is not implemented" );
194
195    /*  newGLCM->numMatrices *= 2;
196        newGLCM->matrixSideLength = maxNumGreyLevels8u*2;
197
198        icvCreateGLCM_Histogram_8uC1R( srcImageStep, srcImageSize, srcImageData,
199                                       newGLCM, numStepDirections,
200                                       stepDirections, memorySteps );
201    */
202    }
203
204    __END__;
205
206    cvFree( &memorySteps );
207    cvFree( &stepDirections );
208
209    if( cvGetErrStatus() < 0 )
210    {
211        cvFree( &newGLCM );
212    }
213
214    return newGLCM;
215}
216
217
218CV_IMPL void
219cvReleaseGLCM( CvGLCM** GLCM, int flag )
220{
221    CV_FUNCNAME( "cvReleaseGLCM" );
222
223    __BEGIN__;
224
225    int matrixLoop;
226
227    if( !GLCM )
228        CV_ERROR( CV_StsNullPtr, "" );
229
230    if( *GLCM )
231        EXIT; // repeated deallocation: just skip it.
232
233    if( (flag == CV_GLCM_GLCM || flag == CV_GLCM_ALL) && (*GLCM)->matrices )
234    {
235        for( matrixLoop = 0; matrixLoop < (*GLCM)->numMatrices; matrixLoop++ )
236        {
237            if( (*GLCM)->matrices[ matrixLoop ] )
238            {
239                cvFree( (*GLCM)->matrices[matrixLoop] );
240                cvFree( (*GLCM)->matrices + matrixLoop );
241            }
242        }
243
244        cvFree( &((*GLCM)->matrices) );
245    }
246
247    if( (flag == CV_GLCM_DESC || flag == CV_GLCM_ALL) && (*GLCM)->descriptors )
248    {
249        for( matrixLoop = 0; matrixLoop < (*GLCM)->numMatrices; matrixLoop++ )
250        {
251            cvFree( (*GLCM)->descriptors + matrixLoop );
252        }
253        cvFree( &((*GLCM)->descriptors) );
254    }
255
256    if( flag == CV_GLCM_ALL )
257    {
258        cvFree( GLCM );
259    }
260
261    __END__;
262}
263
264
265static void
266icvCreateGLCM_LookupTable_8u_C1R( const uchar* srcImageData,
267                                  int srcImageStep,
268                                  CvSize srcImageSize,
269                                  CvGLCM* destGLCM,
270                                  int* steps,
271                                  int numSteps,
272                                  int* memorySteps )
273{
274    int* stepIncrementsCounter = 0;
275
276    CV_FUNCNAME( "icvCreateGLCM_LookupTable_8u_C1R" );
277
278    __BEGIN__;
279
280    int matrixSideLength = destGLCM->matrixSideLength;
281    int stepLoop, sideLoop1, sideLoop2;
282    int colLoop, rowLoop, lineOffset = 0;
283    double*** matrices = 0;
284
285    // allocate memory to the matrices
286    CV_CALL( destGLCM->matrices = (double***)cvAlloc( sizeof(matrices[0])*numSteps ));
287    matrices = destGLCM->matrices;
288
289    for( stepLoop=0; stepLoop<numSteps; stepLoop++ )
290    {
291        CV_CALL( matrices[stepLoop] = (double**)cvAlloc( sizeof(matrices[0])*matrixSideLength ));
292        CV_CALL( matrices[stepLoop][0] = (double*)cvAlloc( sizeof(matrices[0][0])*
293                                                  matrixSideLength*matrixSideLength ));
294
295        memset( matrices[stepLoop][0], 0, matrixSideLength*matrixSideLength*
296                                          sizeof(matrices[0][0]) );
297
298        for( sideLoop1 = 1; sideLoop1 < matrixSideLength; sideLoop1++ )
299        {
300            matrices[stepLoop][sideLoop1] = matrices[stepLoop][sideLoop1-1] + matrixSideLength;
301        }
302    }
303
304    CV_CALL( stepIncrementsCounter = (int*)cvAlloc( numSteps*sizeof(stepIncrementsCounter[0])));
305    memset( stepIncrementsCounter, 0, numSteps*sizeof(stepIncrementsCounter[0]) );
306
307    // generate GLCM for each step
308    for( rowLoop=0; rowLoop<srcImageSize.height; rowLoop++, lineOffset+=srcImageStep )
309    {
310        for( colLoop=0; colLoop<srcImageSize.width; colLoop++ )
311        {
312            int pixelValue1 = destGLCM->forwardLookupTable[srcImageData[lineOffset + colLoop]];
313
314            for( stepLoop=0; stepLoop<numSteps; stepLoop++ )
315            {
316                int col2, row2;
317                row2 = rowLoop + steps[stepLoop*2 + 0];
318                col2 = colLoop + steps[stepLoop*2 + 1];
319
320                if( col2>=0 && row2>=0 && col2<srcImageSize.width && row2<srcImageSize.height )
321                {
322                    int memoryStep = memorySteps[ stepLoop ];
323                    int pixelValue2 = destGLCM->forwardLookupTable[ srcImageData[ lineOffset + colLoop + memoryStep ] ];
324
325                    // maintain symmetry
326                    matrices[stepLoop][pixelValue1][pixelValue2] ++;
327                    matrices[stepLoop][pixelValue2][pixelValue1] ++;
328
329                    // incremenet counter of total number of increments
330                    stepIncrementsCounter[stepLoop] += 2;
331                }
332            }
333        }
334    }
335
336    // normalize matrices. each element is a probability of gray value i,j adjacency in direction/magnitude k
337    for( sideLoop1=0; sideLoop1<matrixSideLength; sideLoop1++ )
338    {
339        for( sideLoop2=0; sideLoop2<matrixSideLength; sideLoop2++ )
340        {
341            for( stepLoop=0; stepLoop<numSteps; stepLoop++ )
342            {
343                matrices[stepLoop][sideLoop1][sideLoop2] /= double(stepIncrementsCounter[stepLoop]);
344            }
345        }
346    }
347
348    destGLCM->matrices = matrices;
349
350    __END__;
351
352    cvFree( &stepIncrementsCounter );
353
354    if( cvGetErrStatus() < 0 )
355        cvReleaseGLCM( &destGLCM, CV_GLCM_GLCM );
356}
357
358
359CV_IMPL void
360cvCreateGLCMDescriptors( CvGLCM* destGLCM, int descriptorOptimizationType )
361{
362    CV_FUNCNAME( "cvCreateGLCMDescriptors" );
363
364    __BEGIN__;
365
366    int matrixLoop;
367
368    if( !destGLCM )
369        CV_ERROR( CV_StsNullPtr, "" );
370
371    if( !(destGLCM->matrices) )
372        CV_ERROR( CV_StsNullPtr, "Matrices are not allocated" );
373
374    CV_CALL( cvReleaseGLCM( &destGLCM, CV_GLCM_DESC ));
375
376    if( destGLCM->optimizationType != CV_GLCM_OPTIMIZATION_HISTOGRAM )
377    {
378        destGLCM->descriptorOptimizationType = destGLCM->numDescriptors = descriptorOptimizationType;
379    }
380    else
381    {
382        CV_ERROR( CV_StsBadFlag, "Histogram-based method is not implemented" );
383//      destGLCM->descriptorOptimizationType = destGLCM->numDescriptors = CV_GLCMDESC_OPTIMIZATION_HISTOGRAM;
384    }
385
386    CV_CALL( destGLCM->descriptors = (double**)
387            cvAlloc( destGLCM->numMatrices*sizeof(destGLCM->descriptors[0])));
388
389    for( matrixLoop = 0; matrixLoop < destGLCM->numMatrices; matrixLoop ++ )
390    {
391        CV_CALL( destGLCM->descriptors[ matrixLoop ] =
392                (double*)cvAlloc( destGLCM->numDescriptors*sizeof(destGLCM->descriptors[0][0])));
393        memset( destGLCM->descriptors[matrixLoop], 0, destGLCM->numDescriptors*sizeof(double) );
394
395        switch( destGLCM->descriptorOptimizationType )
396        {
397            case CV_GLCMDESC_OPTIMIZATION_ALLOWDOUBLENEST:
398                icvCreateGLCMDescriptors_AllowDoubleNest( destGLCM, matrixLoop );
399                break;
400            default:
401                CV_ERROR( CV_StsBadFlag,
402                "descriptorOptimizationType different from CV_GLCMDESC_OPTIMIZATION_ALLOWDOUBLENEST\n"
403                "is not supported" );
404            /*
405            case CV_GLCMDESC_OPTIMIZATION_ALLOWTRIPLENEST:
406                icvCreateGLCMDescriptors_AllowTripleNest( destGLCM, matrixLoop );
407                break;
408            case CV_GLCMDESC_OPTIMIZATION_HISTOGRAM:
409                if(matrixLoop < destGLCM->numMatrices>>1)
410                    icvCreateGLCMDescriptors_Histogram( destGLCM, matrixLoop);
411                    break;
412            */
413        }
414    }
415
416    __END__;
417
418    if( cvGetErrStatus() < 0 )
419        cvReleaseGLCM( &destGLCM, CV_GLCM_DESC );
420}
421
422
423static void
424icvCreateGLCMDescriptors_AllowDoubleNest( CvGLCM* destGLCM, int matrixIndex )
425{
426    int sideLoop1, sideLoop2;
427    int matrixSideLength = destGLCM->matrixSideLength;
428
429    double** matrix = destGLCM->matrices[ matrixIndex ];
430    double* descriptors = destGLCM->descriptors[ matrixIndex ];
431
432    double* marginalProbability =
433        (double*)cvAlloc( matrixSideLength * sizeof(marginalProbability[0]));
434    memset( marginalProbability, 0, matrixSideLength * sizeof(double) );
435
436    double maximumProbability = 0;
437    double marginalProbabilityEntropy = 0;
438    double correlationMean = 0, correlationStdDeviation = 0, correlationProductTerm = 0;
439
440    for( sideLoop1=0; sideLoop1<matrixSideLength; sideLoop1++ )
441    {
442        int actualSideLoop1 = destGLCM->reverseLookupTable[ sideLoop1 ];
443
444        for( sideLoop2=0; sideLoop2<matrixSideLength; sideLoop2++ )
445        {
446            double entryValue = matrix[ sideLoop1 ][ sideLoop2 ];
447
448            int actualSideLoop2 = destGLCM->reverseLookupTable[ sideLoop2 ];
449            int sideLoopDifference = actualSideLoop1 - actualSideLoop2;
450            int sideLoopDifferenceSquared = sideLoopDifference*sideLoopDifference;
451
452            marginalProbability[ sideLoop1 ] += entryValue;
453            correlationMean += actualSideLoop1*entryValue;
454
455            maximumProbability = MAX( maximumProbability, entryValue );
456
457            if( actualSideLoop2 > actualSideLoop1 )
458            {
459                descriptors[ CV_GLCMDESC_CONTRAST ] += sideLoopDifferenceSquared * entryValue;
460            }
461
462            descriptors[ CV_GLCMDESC_HOMOGENITY ] += entryValue / ( 1.0 + sideLoopDifferenceSquared );
463
464            if( entryValue > 0 )
465            {
466                descriptors[ CV_GLCMDESC_ENTROPY ] += entryValue * log( entryValue );
467            }
468
469            descriptors[ CV_GLCMDESC_ENERGY ] += entryValue*entryValue;
470        }
471
472        if( marginalProbability>0 )
473            marginalProbabilityEntropy += marginalProbability[ actualSideLoop1 ]*log(marginalProbability[ actualSideLoop1 ]);
474    }
475
476    marginalProbabilityEntropy = -marginalProbabilityEntropy;
477
478    descriptors[ CV_GLCMDESC_CONTRAST ] += descriptors[ CV_GLCMDESC_CONTRAST ];
479    descriptors[ CV_GLCMDESC_ENTROPY ] = -descriptors[ CV_GLCMDESC_ENTROPY ];
480    descriptors[ CV_GLCMDESC_MAXIMUMPROBABILITY ] = maximumProbability;
481
482    double HXY = 0, HXY1 = 0, HXY2 = 0;
483
484    HXY = descriptors[ CV_GLCMDESC_ENTROPY ];
485
486    for( sideLoop1=0; sideLoop1<matrixSideLength; sideLoop1++ )
487    {
488        double sideEntryValueSum = 0;
489        int actualSideLoop1 = destGLCM->reverseLookupTable[ sideLoop1 ];
490
491        for( sideLoop2=0; sideLoop2<matrixSideLength; sideLoop2++ )
492        {
493            double entryValue = matrix[ sideLoop1 ][ sideLoop2 ];
494
495            sideEntryValueSum += entryValue;
496
497            int actualSideLoop2 = destGLCM->reverseLookupTable[ sideLoop2 ];
498
499            correlationProductTerm += (actualSideLoop1 - correlationMean) * (actualSideLoop2 - correlationMean) * entryValue;
500
501            double clusterTerm = actualSideLoop1 + actualSideLoop2 - correlationMean - correlationMean;
502
503            descriptors[ CV_GLCMDESC_CLUSTERTENDENCY ] += clusterTerm * clusterTerm * entryValue;
504            descriptors[ CV_GLCMDESC_CLUSTERSHADE ] += clusterTerm * clusterTerm * clusterTerm * entryValue;
505
506            double HXYValue = marginalProbability[ actualSideLoop1 ] * marginalProbability[ actualSideLoop2 ];
507            if( HXYValue>0 )
508            {
509                double HXYValueLog = log( HXYValue );
510                HXY1 += entryValue * HXYValueLog;
511                HXY2 += HXYValue * HXYValueLog;
512            }
513        }
514
515        correlationStdDeviation += (actualSideLoop1-correlationMean) * (actualSideLoop1-correlationMean) * sideEntryValueSum;
516    }
517
518    HXY1 =- HXY1;
519    HXY2 =- HXY2;
520
521    descriptors[ CV_GLCMDESC_CORRELATIONINFO1 ] = ( HXY - HXY1 ) / ( correlationMean );
522    descriptors[ CV_GLCMDESC_CORRELATIONINFO2 ] = sqrt( 1.0 - exp( -2.0 * (HXY2 - HXY ) ) );
523
524    correlationStdDeviation = sqrt( correlationStdDeviation );
525
526    descriptors[ CV_GLCMDESC_CORRELATION ] = correlationProductTerm / (correlationStdDeviation*correlationStdDeviation );
527
528    delete [] marginalProbability;
529}
530
531
532CV_IMPL double cvGetGLCMDescriptor( CvGLCM* GLCM, int step, int descriptor )
533{
534    double value = DBL_MAX;
535
536    CV_FUNCNAME( "cvGetGLCMDescriptor" );
537
538    __BEGIN__;
539
540    if( !GLCM )
541        CV_ERROR( CV_StsNullPtr, "" );
542
543    if( !(GLCM->descriptors) )
544        CV_ERROR( CV_StsNullPtr, "" );
545
546    if( (unsigned)step >= (unsigned)(GLCM->numMatrices))
547        CV_ERROR( CV_StsOutOfRange, "step is not in 0 .. GLCM->numMatrices - 1" );
548
549    if( (unsigned)descriptor >= (unsigned)(GLCM->numDescriptors))
550        CV_ERROR( CV_StsOutOfRange, "descriptor is not in 0 .. GLCM->numDescriptors - 1" );
551
552    value = GLCM->descriptors[step][descriptor];
553
554    __END__;
555
556    return value;
557}
558
559
560CV_IMPL void
561cvGetGLCMDescriptorStatistics( CvGLCM* GLCM, int descriptor,
562                               double* _average, double* _standardDeviation )
563{
564    CV_FUNCNAME( "cvGetGLCMDescriptorStatistics" );
565
566    if( _average )
567        *_average = DBL_MAX;
568
569    if( _standardDeviation )
570        *_standardDeviation = DBL_MAX;
571
572    __BEGIN__;
573
574    int matrixLoop, numMatrices;
575    double average = 0, squareSum = 0;
576
577    if( !GLCM )
578        CV_ERROR( CV_StsNullPtr, "" );
579
580    if( !(GLCM->descriptors))
581        CV_ERROR( CV_StsNullPtr, "Descriptors are not calculated" );
582
583    if( (unsigned)descriptor >= (unsigned)(GLCM->numDescriptors) )
584        CV_ERROR( CV_StsOutOfRange, "Descriptor index is out of range" );
585
586    numMatrices = GLCM->numMatrices;
587
588    for( matrixLoop = 0; matrixLoop < numMatrices; matrixLoop++ )
589    {
590        double temp = GLCM->descriptors[ matrixLoop ][ descriptor ];
591        average += temp;
592        squareSum += temp*temp;
593    }
594
595    average /= numMatrices;
596
597    if( _average )
598        *_average = average;
599
600    if( _standardDeviation )
601        *_standardDeviation = sqrt( (squareSum - average*average*numMatrices)/(numMatrices-1));
602
603    __END__;
604}
605
606
607CV_IMPL IplImage*
608cvCreateGLCMImage( CvGLCM* GLCM, int step )
609{
610    IplImage* dest = 0;
611
612    CV_FUNCNAME( "cvCreateGLCMImage" );
613
614    __BEGIN__;
615
616    float* destData;
617    int sideLoop1, sideLoop2;
618
619    if( !GLCM )
620        CV_ERROR( CV_StsNullPtr, "" );
621
622    if( !(GLCM->matrices) )
623        CV_ERROR( CV_StsNullPtr, "Matrices are not allocated" );
624
625    if( (unsigned)step >= (unsigned)(GLCM->numMatrices) )
626        CV_ERROR( CV_StsOutOfRange, "The step index is out of range" );
627
628    dest = cvCreateImage( cvSize( GLCM->matrixSideLength, GLCM->matrixSideLength ), IPL_DEPTH_32F, 1 );
629    destData = (float*)(dest->imageData);
630
631    for( sideLoop1 = 0; sideLoop1 < GLCM->matrixSideLength;
632                        sideLoop1++, (float*&)destData += dest->widthStep )
633    {
634        for( sideLoop2=0; sideLoop2 < GLCM->matrixSideLength; sideLoop2++ )
635        {
636            double matrixValue = GLCM->matrices[step][sideLoop1][sideLoop2];
637            destData[ sideLoop2 ] = (float)matrixValue;
638        }
639    }
640
641    __END__;
642
643    if( cvGetErrStatus() < 0 )
644        cvReleaseImage( &dest );
645
646    return dest;
647}
648
649