Mercurial > hg > ltpda

comparison src/ltpda_polyreg/ltpda_polyreg.c @ 0:f0afece42f48

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Import.
author Daniele Nicolodi <nicolodi@science.unitn.it>
date Wed, 23 Nov 2011 19:22:13 +0100
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:f0afece42f48
1 /*
2 Mex file that implements polynomial detrending of a data vector.
3
4 M Hewitson 5-02-08
5
6 $Id: ltpda_polyreg.c,v 1.4 2009/12/17 08:02:08 mauro Exp $
7 */
8
9
10 #include <math.h>
11 #include <stdio.h>
12 #include <stdlib.h>
13 #include <string.h>
14 #include <mex.h>
15
16 #include "ltpda_polyreg.h"
17 #include "version.h"
18 #include "../c_sources/polyreg.c"
19
20 #define DEBUG 0
21
22 /*
23 Matlab mex file to make polynomial detrending of a data vector
24
25 * Inputs
26 * - data
27 * - N
28 *
29 function [y,a] = ltpda_polyreg(x, N);
30
31 */
32 void mexFunction( int nlhs, mxArray *plhs[],
33 int nrhs, const mxArray *prhs[])
34 {
35 /* Parse inputs */
36 if( (nrhs == 0) && (nlhs == 0) )
37 {
38 print_usage(VERSION);
39 return;
40 }
41 else if ( (nrhs == 2) && (nlhs >= 1) ) /* let's go */
42 {
43 double *xdata;
44 long int nData;
45 int order;
46 double *yptr, *aptr;
47
48
49 /* Extract inputs */
50 xdata = mxGetPr(prhs[0]); /* Pointer to data */
51 nData = mxGetNumberOfElements(prhs[0]); /* Number of data points */
52 order = (int)mxGetScalar(prhs[1]); /* Order of detrending */
53
54 if (order > 10 || order < -1)
55 mexErrMsgTxt("Detrending order must be between -1 and 10");
56
57 /* Set output matrices */
58 plhs[0] = mxCreateDoubleMatrix(nData, 1, mxREAL);
59
60 /* Get pointers to output matrices */
61 yptr = mxGetPr(plhs[0]);
62 if (nlhs == 2)
63 {
64 plhs[1] = mxCreateDoubleMatrix(order+1, 1, mxREAL);
65 aptr = mxGetPr(plhs[1]);
66 }
67 else
68 {
69 aptr = (double*)mxCalloc(order+1, sizeof(double)); /* detrending output */
70 }
71
72 /* Detrend segment */
73 switch (order)
74 {
75 case -1:
76 /* no detrending */
77 memcpy(yptr, xdata, nData*sizeof(double));
78 break;
79 case 0:
80 /* mean removal */
81 polyreg0 (xdata, nData, yptr, aptr);
82 break;
83 case 1:
84 /* linear detrending */
85 polyreg1 (xdata, nData, yptr, aptr);
86 break;
87 case 2:
88 /* 2nd order detrending */
89 polyreg2 (xdata, nData, yptr, aptr);
90 break;
91 case 3:
92 /* 3rd order detrending */
93 polyreg3 (xdata, nData, yptr, aptr);
94 break;
95 case 4:
96 /* 4th order detrending */
97 polyreg4 (xdata, nData, yptr, aptr);
98 break;
99 case 5:
100 /* 5th order detrending */
101 polyreg5 (xdata, nData, yptr, aptr);
102 break;
103 case 6:
104 /* 6th order detrending */
105 polyreg6 (xdata, nData, yptr, aptr);
106 break;
107 case 7:
108 /* 7th order detrending */
109 polyreg7 (xdata, nData, yptr, aptr);
110 break;
111 case 8:
112 /* 8th order detrending */
113 polyreg8 (xdata, nData, yptr, aptr);
114 break;
115 case 9:
116 /* 9th order detrending */
117 polyreg9 (xdata, nData, yptr, aptr);
118 break;
119 case 10:
120 /* 10th order detrending */
121 polyreg10(xdata, nData, yptr, aptr);
122 break;
123 }
124
125 /*// if (nlhs == 2)
126 // {
127 // // Free coefficient array
128 // mxFree(aptr);
129 // }*/
130 }
131 else /* we have an error */
132 {
133 print_usage(VERSION);
134 mexErrMsgTxt("### incorrect usage");
135 }
136 }
137
138 /*
139 * Output usage to MATLAB terminal
140 *
141 */
142 void print_usage(char *version)
143 {
144 mexPrintf("ltpda_polyreg version %s\n", version);
145 mexPrintf(" usage: function [y, a] = ltpda_polyreg(x, order); \n");
146 }
147
148