fractional_diff.c

来自「最大似然估计算法」· C语言 代码 · 共 46 行

#include "timeseries.h"

void fractional_diff(double *a, int n, double degree) {

	/************************************************************************/
	/*    fractional_diff(double *a, int n, double degree)                  */
	/*                                                                      */
	/*    Generate a differencing/integrating matrix of any order including */
	/*    fractional.  These are the discrete-time versions of fractional   */
	/*    differentiation and integration.                                  */
	/*                                                                      */
	/*    length    size of filter matrix                                   */
	/*    degree    differencing/integrating exponent:                      */
	/*              differencing => degree > 0                              */
	/*              integrating  => degree < 0                              */
	/*                                                                      */
	/*    matrix    square filter matrix                                    */
	/*                                                                      */
	/* See Hosking, R. M., 1981, Biometrika, 68, pp. 165-176                */
	/*                                                                      */
	/* Hadley Johnson, [e-mail], 28-May-1998 version               */
	/* Simon Williams, [e-mail], 08-Nov-2000 c-version               */
	/************************************************************************/

	int i, j;

	double *b, *c;

	b = (double *)calloc((size_t)n, sizeof(double));
	c = (double *)calloc((size_t)n, sizeof(double));

	b[0] = 1.0;
	
	for (j = 0; j < n-1; j++)  b[j+1] = ((double) j - degree) / ((double) (j+1));

	c[0] = b[0];
	for (j = 1; j < n; j++) c[j] = c[j-1] * b[j];
	
	for (j = 0; j < n; j++) {
		for (i = j; i < n; i++) a[i + j * n] = c[i-j];
	}

	free(b);
	free(c);
}