Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Perform autocorrelation with vDSP_conv from Apple Accelerate Framework

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

like image 937
MrHappyAsthma Avatar asked Jun 06 '12 16:06

MrHappyAsthma


1 Answers

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);
like image 199
Rasman Avatar answered Nov 19 '22 03:11

Rasman