Perform autocorrelation with vDSP_conv from Apple Accelerate Framework
Asked Answered
L

1

7

I need to perform the autocorrelation of an array (vector) but I am having trouble finding the correct way to do so. I believe that I need the method "vDSP_conv" from the Accelerate Framework, but I can't follow how to successfully set it up. The thing throwing me off the most is the need for 2 inputs. Perhaps I have the wrong function, but I couldn't find one that operated on a single vector.

The documentation can be found here

Copied from the site

vDSP_conv

Performs either correlation or convolution on two vectors; single precision.

void vDSP_conv ( const float __vDSP_signal[], vDSP_Stride __vDSP_signalStride, const float __vDSP_filter[], vDSP_Stride __vDSP_strideFilter, float __vDSP_result[], vDSP_Stride __vDSP_strideResult, vDSP_Length __vDSP_lenResult, vDSP_Length __vDSP_lenFilter );

Parameters

__vDSP_signal

Input vector A. The length of this vector must be at least __vDSP_lenResult + __vDSP_lenFilter - 1.

__vDSP_signalStride

The stride through __vDSP_signal.

__vDSP_filter

Input vector B.

__vDSP_strideFilter

The stride through __vDSP_filter.

__vDSP_result

Output vector C.

__vDSP_strideResult

The stride through __vDSP_result.

__vDSP_lenResult

The length of __vDSP_result.

__vDSP_lenFilter

The length of __vDSP_filter.

For an example, just assume you have an array of float x = [1.0, 2.0, 3.0, 4.0, 5.0]. How would I take the autocorrelation of that?

The output should be something similar to float y = [5.0, 14.0, 26.0, 40.0, 55.0, 40.0, 26.0, 14.0, 5.0] //generated using Matlab's xcorr(x) function

Longitudinal answered 6/6, 2012 at 16:7 Comment(0)
V
4

performing autocorrelation simply means you take the cross-correlation of one vector with itself. There is nothing fancy about it.

so in your case, do:

vDSP_conv(x, 1, x, 1, result, 1, 2*len_X-1, len_X); 

check a sample code for more details: (which does a convolution)

http://disanji.net/iOS_Doc/#documentation/Performance/Conceptual/vDSP_Programming_Guide/SampleCode/SampleCode.html

EDIT: This borders on ridiculous, but you need to offset the x value by a specific number of zeros, which is just crazy.

the following is a working code, just set filter to the value of x you desire, and it will put the rest in the correct position:

float          *signal, *filter, *result;

int32_t         signalStride, filterStride, resultStride;

uint32_t        lenSignal, filterLength, resultLength;

uint32_t        i;



filterLength = 5;

resultLength = filterLength*2 -1;

lenSignal = ((filterLength + 3) & 0xFFFFFFFC) + resultLength;



signalStride = filterStride = resultStride = 1;



printf("\nConvolution ( resultLength = %d, "

       "filterLength = %d )\n\n", resultLength, filterLength);



/* Allocate memory for the input operands and check its availability. */

signal = (float *) malloc(lenSignal * sizeof(float));

filter = (float *) malloc(filterLength * sizeof(float));

result = (float *) malloc(resultLength * sizeof(float));



for (i = 0; i < filterLength; i++)

    filter[i] = (float)(i+1);

for (i = 0; i < resultLength; i++)
    if (i >=resultLength- filterLength)
        signal[i] = filter[i - filterLength+1];


/* Correlation. */

vDSP_conv(signal, signalStride, filter, filterStride,

          result, resultStride, resultLength, filterLength);


printf("signal: ");
for (i = 0; i < lenSignal; i++)        
    printf("%2.1f ", signal[i]);


printf("\n filter: ");
for (i = 0; i < filterLength; i++)
    printf("%2.1f ", filter[i]);

printf("\n result: ");
for (i = 0; i < resultLength; i++)
    printf("%2.1f ", result[i]);


/* Free allocated memory. */

free(signal);

free(filter);

free(result);
Voiceful answered 6/6, 2012 at 18:31 Comment(6)
I just tried this with my example. But it outputs slightly wrong data. What am I doing wrong? x = 1 to 5 and result is size 9 (5 * 2 -1) vDSP_conv(x, 1, x, 1, result, 1, 9, 5); -> outputs -> 55.0, 40.0, 26.0, 14.0, 5.0, 0.0, 0.0, 0.0, 294.0Longitudinal
I think I made a mistake, you may need to pad x with zeroes... I don't have my Mac in front of me, but it seems that the input should be longer then the output (which is strange to say the least)Voiceful
Blah :P That stinks. Ill try messing around with it. If you come across a way to fix it for sure then please post :DLongitudinal
@Longitudinal it took awhile but I figured it out. See above edit. Good luckVoiceful
That is so weird that it has to be done like that, but it does work! :D Thank you so much!!!Longitudinal
This is an old question - but in case it is of use - if you want the linear autocorrelation then you need to pad the input with zeros of the same length as the data - if you don't then you'll get a cyclic autocorrelation and hence slightly wrong answers - I guess because internally this is using FFTs.Fourteenth

© 2022 - 2024 — McMap. All rights reserved.