1 |
// $Header: $ |
2 |
|
3 |
#include <stdio.h> |
4 |
#include <stdlib.h> |
5 |
#include <assert.h> |
6 |
#include <string.h> |
7 |
#include <complex.h> |
8 |
#include <math.h> |
9 |
#include "ctypes.h" |
10 |
#include "detrend_C99.h" |
11 |
//#include "xmem.h" |
12 |
|
13 |
//#define DEBUGFILE //uncomment this to write debugging files |
14 |
|
15 |
#ifdef FLOAT |
16 |
#define TYPE FLOAT |
17 |
#include "detrend_code_C99.h" |
18 |
#undef TYPE |
19 |
#endif |
20 |
|
21 |
#ifdef DOUBLE |
22 |
#define TYPE DOUBLE |
23 |
#include "detrend_code_C99.h" |
24 |
#undef TYPE |
25 |
#endif |
26 |
|
27 |
#ifdef COMPLEXFLOAT |
28 |
static int count2=0; |
29 |
void cdetrend_discontig( int n, _Complex float *data, int *isgood, |
30 |
int degree, int length, int skip, |
31 |
int m, int *sect, int detrend_first) |
32 |
{ |
33 |
int i,first,last; |
34 |
// Seperately detrend m sections of data separated by |
35 |
// discontinuities. |
36 |
// The first section begins on data[0], |
37 |
// the m-1 first sections end on data[sect_last[i]], i=0,1,...,m-2, |
38 |
// and the last section ends on data[n-1]. |
39 |
|
40 |
if (m==0) |
41 |
cdetrend(n,data,isgood,degree,length,skip,detrend_first); |
42 |
else |
43 |
{ |
44 |
for (i=0; i<m; i++) |
45 |
{ |
46 |
first = sect[2*i]; |
47 |
last = sect[2*i+1]; |
48 |
#ifdef DEBUGOUT |
49 |
printf("Detrending [%d:%d]\n",first,last); |
50 |
#endif |
51 |
cdetrend(last-first+1,&data[first],&isgood[first],degree,length,skip,detrend_first); |
52 |
} |
53 |
} |
54 |
} |
55 |
|
56 |
void cdetrend( int n, _Complex float *data, int *isgood, int degree, |
57 |
int length, int skip, int detrend_first) |
58 |
{ |
59 |
int i; |
60 |
|
61 |
float *real, *imag; |
62 |
real = (float *) malloc(n*sizeof(float)); |
63 |
for (i=0;i<n;i++) |
64 |
real[i] = crealf(data[i]); |
65 |
sdetrend( n, real, isgood, degree, length, skip, detrend_first); |
66 |
|
67 |
imag = (float *) malloc(n*sizeof(float)); |
68 |
for (i=0;i<n;i++) |
69 |
imag[i] = cimagf(data[i]); |
70 |
sdetrend( n, imag, isgood, degree, length, skip, detrend_first); |
71 |
for (i=0;i<n;i++) |
72 |
data[i] = real[i] + _Complex_I*imag[i]; |
73 |
|
74 |
#ifdef DEBUGFILE |
75 |
{ |
76 |
FILE *fh; |
77 |
char name[1024]; |
78 |
printf("writing debug2 file #%d\n",count2); |
79 |
sprintf(name,"detrend_debug2_%d.out",count2++); |
80 |
fh = fopen(name,"w"); |
81 |
assert(fh); |
82 |
for (i=0; i<n; i++) |
83 |
fprintf(fh,"%e %e\n",crealf(data[i]),cimagf(data[i])); |
84 |
fclose(fh); |
85 |
} |
86 |
#endif |
87 |
|
88 |
free(real); |
89 |
free(imag); |
90 |
} |
91 |
#endif |
92 |
|
93 |
|
94 |
#ifdef COMPLEXDOUBLE |
95 |
void zdetrend_discontig( int n, _Complex double *data, int *isgood, |
96 |
int degree, int length, int skip, |
97 |
int m, int *sect, int detrend_first) |
98 |
{ |
99 |
int i,first,last; |
100 |
// Seperately detrend m sections of data separated by |
101 |
// discontinuities. |
102 |
// The first section begins on data[0], |
103 |
// the m-1 first sections end on data[sect_last[i]], i=0,1,...,m-2, |
104 |
// and the last section ends on data[n-1]. |
105 |
|
106 |
if (m==0) |
107 |
zdetrend(n,data,isgood,degree,length,skip,detrend_first); |
108 |
else |
109 |
{ |
110 |
first = 0; |
111 |
for (i=0; i<m; i++) |
112 |
{ |
113 |
first = sect[2*i]; |
114 |
last = sect[2*i+1]; |
115 |
zdetrend(last-first+1,&data[first],&isgood[first],degree,length,skip,detrend_first); |
116 |
} |
117 |
} |
118 |
} |
119 |
|
120 |
|
121 |
void zdetrend( int n, _Complex double *data, int *isgood, int degree, |
122 |
int length, int skip, int detrend_first) |
123 |
{ |
124 |
int i; |
125 |
double *real, *imag; |
126 |
|
127 |
real = (double *) malloc(n*sizeof(double)); |
128 |
for (i=0;i<n;i++) |
129 |
real[i] = creal(data[i]); |
130 |
ddetrend( n, real, isgood, degree, length, skip, detrend_first); |
131 |
|
132 |
imag = (double *) malloc(n*sizeof(double)); |
133 |
for (i=0;i<n;i++) |
134 |
imag[i] = cimag(data[i]); |
135 |
ddetrend( n, imag, isgood, degree, length, skip, detrend_first); |
136 |
for (i=0;i<n;i++) |
137 |
data[i] = real[i] + _Complex_I*imag[i]; |
138 |
|
139 |
free(real); |
140 |
free(imag); |
141 |
} |
142 |
#endif |