转载自微信公众号
原创 逸珺 嵌入式客栈
#include <stdio.h>
#include <math.h>
/* 返回值在区间: [-1,1] */
/* 如返回-10,则证明输入参数无效 */
#define delta 0.0001f
double calculate_corss_correlation(double *s1, double *s2,int n)
{
double sum_s12 = 0.0;
double sum_s1 = 0.0;
double sum_s2 = 0.0;
double sum_s1s1 = 0.0; //s1^2
double sum_s2s2 = 0.0; //s2^2
double pxy = 0.0;
double temp1 = 0.0;
double temp2 = 0.0;
if( s1==NULL || s2==NULL || n<=0)
return -10;
for(int i=0;i<n;i++)
{
sum_s12 += s1[i]*s2[i];
sum_s1 += s1[i];
sum_s2 += s2[i];
sum_s1s1 += s1[i]*s1[i];
sum_s2s2 += s2[i]*s2[i];
}
temp1 = n*sum_s1s1-sum_s1*sum_s1;
temp2 = n*sum_s2s2-sum_s2*sum_s2;
/* 分母不可为0 */
if( (temp1>-delta && temp1<delta) ||
(temp2>-delta && temp2<delta) ||
(temp1*temp2<=0) )
{
return -10;
}
pxy = (n*sum_s12-sum_s1*sum_s2)/sqrt(temp1*temp2);
return pxy;
}
double s1[30] = {
0.309016989,0.587785244,0.809016985,0.95105651,1,0.951056526,
0.809017016,0.587785287,0.30901704,5.35898E-08,0,0,
0,0,0,0,0,0,
0,0,0,0,0,0,
0,0,0,0,0,0
};
double s2[30] = {
0.343282816,0.686491368,0.874624132,0.99459642,1.008448609,
1.014252458,0.884609221,0.677632906,0.378334666,0.077878732,
0.050711886,0.066417083,0.088759401,0.005440732,0.04225661,
0.035349939,0.0631196,0.007566056,0.053183895,0.073143706,
0.080285063,0.030110227,0.044781145,0.01875573,0.08373928,
0.04550342,0.038880858,0.040611891,0.046116826,0.087670453
};
int main(void)
{
double pxy;
double s3[30];
pxy = calculate_corss_correlation(s1,s2,30);
printf("pxy of s1 and s2:%f\n",pxy);
pxy = calculate_corss_correlation(s1,s1,30);
printf("pxy of s1 and s1:%f\n",pxy);
for(int i=0;i<n;i++)
{
s3[i] = -1*s1[i];
}
pxy = calculate_corss_correlation(s1,s3,30);
printf("pxy of s1 and s3:%f\n",pxy);
return 0;
}
pxy of s1 and s2:0.997435
pxy of s1 and s1:1.000000
pxy of s1 and s1:-1.000000