Commit | Line | Data |
---|---|---|
e161499b BA |
1 | #include <stdlib.h> |
2 | #include <R.h> | |
3 | #include <Rinternals.h> | |
4 | #include <Rmath.h> | |
5 | ||
6 | #include <stdio.h> | |
7 | ||
8 | SEXP filter(SEXP cwt_) | |
9 | { | |
10 | int L = INTEGER(getAttrib(cwt_, R_DimSymbol))[0], | |
11 | D = INTEGER(getAttrib(cwt_, R_DimSymbol))[1]; | |
12 | double *cwt = REAL(cwt_); | |
13 | SEXP R_fcwt; | |
14 | PROTECT(R_fcwt = allocMatrix(REALSXP, L, D)); | |
15 | double* fcwt = REAL(R_fcwt); | |
16 | ||
17 | //TODO: coding style is terrible... no time for now. | |
18 | for (int col=0; col<D; col++) | |
19 | { | |
20 | double v1 = cwt[0]; | |
21 | double ma = v1 + cwt[1] + cwt[2]; | |
22 | for (int i=1; i<L-2; i++) | |
23 | { | |
24 | fcwt[i] = ma / 3.; | |
25 | ma = ma - v1 + cwt[i+2]; | |
26 | v1 = cwt[i]; | |
27 | } | |
28 | fcwt[0] = cwt[0]; | |
29 | fcwt[L-2] = ma / 3.; | |
30 | fcwt[L-1] = cwt[L-1]; | |
31 | cwt = cwt + L; | |
32 | fcwt = fcwt + L; | |
33 | } | |
34 | ||
35 | UNPROTECT(1); | |
36 | return R_fcwt; | |
37 | } |