<commit_before>/*
 *  RecursiveFilter.cpp
 *
 * Copyright 2011 <NAME>
 *  All rights reserved.
 *
 */

#include "RecursiveFilter.h"
#include <iostream>
#include <cmath>

RecursiveFilter::RecursiveFilter(const int& fOrder, const int& fLengthScale)
{
	
	// Create a filter based on the specified order and lengthscale
	order = fOrder;
	lengthScale = fLengthScale;
	getFilterCoefficients();
	
}

RecursiveFilter::~RecursiveFilter()
{
}

void RecursiveFilter::getFilterCoefficients()
{
	real sigma = (lengthScale*lengthScale)/2;
	real b[5][5];
	real a[201][201];
	real p[201];
	real K[9];
	real D[9];
	
	// b coefficients
	b[1][1] = 1;
	b[1][2] = 1/12;
	b[1][3] = 1/90;
	b[1][4] = 1/560;
	b[2][2] = 1;
	b[2][3] = 1/6;
	b[2][4] = 7/240;
	b[3][3] = 1;
	b[3][4] = 1/4;
	b[4][4] = 1;
	
    // Compute polynomial coefficients
	K[2] = b[1][1]*sigma;
	K[4] = b[1][2]*sigma + b[2][2]*(pow(sigma,2))/factorial(2);
	K[6] = b[1][3]*sigma + b[2][3]*(pow(sigma,2))/factorial(2)
		+ b[3][3]*(pow(sigma,3))/factorial(3);
	K[8] = b[1][4]*sigma + b[2][4]*(pow(sigma,2))/factorial(2)
		+ b[3][4]*(pow(sigma,3))/factorial(3) + b[4][4]*(pow(sigma,4))/factorial(4);
	
	for (int i=0;i<=order*2;i++) {
		D[i] = 0;
	}
	
	for (int n=2;n<=order*2;n+=2) {
		for (int i=0;i<=n;i++) {
			int index = n/2 - i;
			real coeff = (pow((-1),(n/2)))*(pow((-1),i))*factorial(n)/(factorial(n-i)*factorial(i));
			D[index+order] += K[n]*coeff;
		}
	}
	D[order] += 1;
	
	for (int i=0;i<=200;i++) {
		for (int j=0;j<=200;j++) {
			if (abs(j-i) <= order) {
				a[i][j] = D[j-i+order];
			} else {
				a[i][j] = 0;
			}
		}
	}
	
	// Cholesky decomp
	for (int i=0;i<=200;i++) {
		for (int j=0;j<=200;j++) {
			real sum=a[i][j];
			for (int k=i-1;k>=0;k--) {
				sum -= a[i][k]*a[j][k];
			}
			if (i == j) {
				if (sum <= 0.0) { 
					std::cout << "cholesky failed at i,j sum\n";
					break;
				} else {
					p[i] = sqrt(sum);
				}
			} else {<commit_msg>Fix abs command to use std::abs<commit_after>/*
 *  RecursiveFilter.cpp
 *
 * Copyright 2011 <NAME>
 *  All rights reserved.
 *
 */

#include "RecursiveFilter.h"
#include <iostream>
#include <cmath>

RecursiveFilter::RecursiveFilter(const int& fOrder, const int& fLengthScale)
{
	
	// Create a filter based on the specified order and lengthscale
	order = fOrder;
	lengthScale = fLengthScale;
	getFilterCoefficients();
	
}

RecursiveFilter::~RecursiveFilter()
{
}

void RecursiveFilter::getFilterCoefficients()
{
	real sigma = (lengthScale*lengthScale)/2;
	real b[5][5];
	real a[201][201];
	real p[201];
	real K[9];
	real D[9];
	
	// b coefficients
	b[1][1] = 1;
	b[1][2] = 1/12;
	b[1][3] = 1/90;
	b[1][4] = 1/560;
	b[2][2] = 1;
	b[2][3] = 1/6;
	b[2][4] = 7/240;
	b[3][3] = 1;
	b[3][4] = 1/4;
	b[4][4] = 1;
	
    // Compute polynomial coefficients
	K[2] = b[1][1]*sigma;
	K[4] = b[1][2]*sigma + b[2][2]*(pow(sigma,2))/factorial(2);
	K[6] = b[1][3]*sigma + b[2][3]*(pow(sigma,2))/factorial(2)
		+ b[3][3]*(pow(sigma,3))/factorial(3);
	K[8] = b[1][4]*sigma + b[2][4]*(pow(sigma,2))/factorial(2)
		+ b[3][4]*(pow(sigma,3))/factorial(3) + b[4][4]*(pow(sigma,4))/factorial(4);
	
	for (int i=0;i<=order*2;i++) {
		D[i] = 0;
	}
	
	for (int n=2;n<=order*2;n+=2) {
		for (int i=0;i<=n;i++) {
			int index = n/2 - i;
			real coeff = (pow((-1),(n/2)))*(pow((-1),i))*factorial(n)/(factorial(n-i)*factorial(i));
			D[index+order] += K[n]*coeff;
		}
	}
	D[order] += 1;
	
	for (int i=0;i<=200;i++) {
		for (int j=0;j<=200;j++) {
			if (std::abs((double)(j-i)) <= order) {
				a[i][j] = D[j-i+order];
			} else {
				a[i][j] = 0;
			}
		}
	}
	
	// Cholesky decomp
	for (int i=0;i<=200;i++) {
		for (int j=0;j<=200;j++) {
			real sum=a[i][j];
			for (int k=i-1;k>=0;k--) {
				sum -= a[i][k]*a[j][k];
			}
			if (i == j) {
				if (sum <= 0.0) { 
					std::cout << "cholesky failed at i,j sum\n";
					break;
				} else {
					p[i] = sqrt(sum);
				}
			} else {
{'fraction_non_alphanumeric': 0.2590189462756294, 'fraction_numerical': 0.06696080975862964, 'mean_word_length': 3.4503464203233256, 'pattern_counts': {'":': 0, '<': 31, '<?xml version=': 0, '>': 11, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '2259531', 'n_tokens_mistral': 2102, 'n_tokens_neox': 1849, 'n_words': 321}